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^ ' The calculation of electronic properties of materials is an important task of solid state theory, 
^^Ibeit particularly difficult if electronic correlations are strong, for example in transition metals, 

I tih eir oxides and in /-electron systems. The standard approach to material calculations, the density 
functional theory in its local density approximation (LDA), incorporates electronic correlations only 

^s^ery rudimentarily and fails if the correlations are strong. Encouraged by the success of dynamical 
K^nean field theory (DMFT) in dealing with strongly correlated model Hamiltonians, physicists from 

^..the bandstructure and the many-body community have joined forces and developed a combined 

■^/LDA-I-DMFT method recently. Depending on the strength of electronic correlations, this new ap- 

^Jproach yields a weakly correlated metal as in LDA, a strongly correlated metal, or a Mott insulator. 

CNBy now, this approach is widely regarded as a breakthrough for electronic structure calculations of 

^"^trongly correlated materials. 

T~H The author will review this LDA+DMFT method and also discuss alternative approaches to 
l/~]pmploy DMFT in electronic structure calculations, for example, by replacing the LDA part by the 
f~^ o-called GW approximation. Different methods to solve the DMFT equations are introduced with 
■"""^^ focus on those that are suitable for realistic calculations with many orbitals. An overview of the 
■^liuccessful application of LDA+DMFT to a wide variety of materials, ranging from Pu and Ce, to Fe 
and Ni, to numerous transition metal oxides, is given. 
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1 Introduction 

One of the most important challenges of theoretical physics is the develop- 
ment of reliable methods for the quantitative calculation of material prop- 
erties. In solid state theory, we know the Hamiltonian to do these material 
calculations. For example, if we neglect relativistic corrections and employ 
the Born-Oppenheimer (1) approximation, this Hamiltonian consists of three 
terms: the kinetic energy, the lattice potential and the Coulomb interaction 
between the electrons: 



2me 4^ 47reo |rj - R; 



1 >p ^ 1 

^ 2 ^ 47reo \ri - rj\' ^ ' 



Here, and R; denote the position of electron i and ion I with charge — e and 
Zie, respectively; A.j is the Laplace operator for the kinetic energy of electrons 
with mass nie', cq and h are the vacuum dielectric and Planck constant. 
The three terms of Hamiltonian [1] describe the movement of the electrons. 
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correlations 



Solid State Hamiltonian: 

To calculate the physical properties 
of a given material, one has to take 
into account three terms (Hamilto- 
nian [1]): The kinetic energy due 
to which the electrons move (ar- 
row), the lattice potential of the 
ions and the Coulomb interaction 
between the electrons. Due to the 
latter, the moving electron is re- 
pelled by the electron already lo- 
cated at this site. It is energetically 
favourable if the depicted electron 
hops somewhere else, as indicated 
by the arrow. Hence, the movement 
of every electron is correlated with 
that of every other. 




time-averaged 
electron density 



attice potential 



LDA approximation: 

LDA is an approximation which 
allows for calculating material 
properties, but dramatically sim- 
plifies the electronic correlations: 
The LDA bandstructure corre- 
sponds to the simplification that 
every electron moves indepen- 
dently, i.e. uncorrelated, within 
a time-averaged local density of 
the other electrons (visualised as 
static clouds in the figure). 



the attractive lattice potential of the (fixed) ions and the mutual Coulomb 
repulsion of the electrons, see Fig. 1 for an illustration. But while we know 
the Hamiltonian, we cannot solve it, even not numerically, if more than a 
very few [0(10)] electrons are involved. This problem is due to the last term, 
the Coulomb interaction, which correlates the movement of every electron i 
with every other electron j. The numerical effort to solve Eq. [1] quantum- 
mechanically grows exponentially with the number of correlated electrons. 
Such an exponential problem cannot be solved for significantly many electrons, 
even if computer power continued to grow rapidly. 

In this situation, we have two possibilities: Either we dramatically simplify 
the Hamiltonian [1], hoping that the simplified model allows for a qualitative 
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understanding including electronic correlation. Or we employ equally dramatic 
approximations to deal with Eq. [1] directly. These two strategies have been 
followed by the two large communities of solid state theory: the many-body 
model Hamiltonian and the density functional community. 

Within the density functional theory (DFT) introduced by Hohenberg and 
Kohn (2), the local density approximation (LDA) turned out to be unex- 
pectedly successful, and established itself as the method for realistic solid 
state calculations in the last century; for reviews see e.g. Refs. (3; 4; 5; 6). 
This is surprising because LDA represents a substantial approximation to the 
Coulomb interaction. Basically, an electron at sees a time-averaged local 
density p{vi) of the other electrons, with a corresponding local LDA potential 
^LDA(p(ri)), as visualised in Fig. 1. This reduces the many-body problem to 
a single electron calculation in a local potential. 

The success of LDA shows that this treatment is actually sufficient for many 
materials, both for calculating ground state energies and bandstructures, im- 
plying that electronic correlations are rather weak in these materials. But, 
there are important classes of materials where LDA fails, such as transition 
metal oxides or heavy Fermion systems. For example, LDA predicts La2Cu04, 
Cr-doped V2O3 and NiO to be metals above the antiferromagnetic ordering 
temperature whereas, in reality, these materials are insulators (7; 8; 9; 10; 11). 
Similarly, the LDA bandstructure is in strong disagreement with experiment 
for /-electron systems. In these two classes of materials, the valence electrons 
partially occupy the d or / orbitals. If there are two electrons in these narrow 
orbitals on the same lattice site, the distance |rj — rj| is particularly short, 
and electronic correlations are thus particularly strong. Hence, the approxi- 
mate LDA treatment of Hamiltonian [1] fails. 

To describe these d- and /-electron systems correctly, genuine many-body 
effects have to be taken into account. The transition metal oxides mentioned 
above are, for example, Mott insulators (12; 13; 14). Since the d orbitals are 
only partially filled the LDA bandstructure predicts metallic behaviour, see 
the left panel of Fig. 2. But instead, the on- (lattice-) site Coulomb repulsion 
U splits the metallic LDA bands into two sets of Hubbard (15) bands, as 
in the right panel of Fig. 2. Let us for a moment assume that the average 
number of d or / electrons is one per lattice site ('integer filling'). Then, one 
can envisage the lower Hubbard band as consisting of all states with exactly 
one electron on every lattice site and the upper Hubbard band as those states 
where two electrons are in d or / orbitals on the same lattice site. Since it costs 
an energy U to have two electrons on the same lattice sites, the latter states 
are completely empty and the former completely filled with a gap of size U in 
between. Such kind of Mott-Hubbard physics and the associated energy gain 
is completely missing in the LDA. 

This shortcoming can be overcome by the self interaction correction (SIC) 
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LDA+DMFT 



LDA+U 



Weakly 

correlated metal: 

The on-site 
Coulomb inter- 
action is weak 
compared to the 
LDA bandwidth: 
U <^W; and LDA 
gives the correct 
answer, typically a 
(weakly correlated) 
metal for which we 
schematically draw 
a density of states 
at energies E and 
the bandstructure, 
i.e. E vs. wave 
vector k. 



U/W 

Figure 2. 

Strongly 

correlated metal: 

In this intermediate 
regime, one has already 
Hubbard bands, like for 
U/W > 1 (right hand 
side), but at the same time 
a remainder of the weakly 
correlated LDA metal (left 
hand side), in form of a 
quasiparticle peak: The 
(U = 0) LDA bandstructure 
is reproduced, albeit with 
its width and weight re- 
duced by a factor Z and life 
time effects which result in 
a Lorentzian broadening of 
the quasiparticle levels. 



Mott insulator: 

If the Coulomb in- 
teraction U becomes 
large (C/ > VF), the 
LDA band splits 
into two Hubbard 
bands, and we have 
a Mott insulator 
with only the lower 
band occupied (at 
integer fillings). Such 
a splitting can be 
described by the 
so-cahed LDA-hU 
method with the 
drawbacks discussed 
in the text. 



The entire parameter regime is described by LDA-I-DMFT. 



to LDA by Perdew and Zunger (16; 17; 18) or by the LDA-I-U method of 
Anisimov et al. (19) which we will discuss in Section 5.1. In the presence of 
orbital or magnetic ordering these approaches yield insulating spectra, similar 
to that of the paramagnetic Mott insulator displayed in the right panel of Fig. 
2. Important differences are that the LDA-I-U spectra are completely coherent 
and that the LDA-I-U solution has not the high entropy of the finite temper- 
ature paramagnetic solution, because of the orbital or magnetic ordering. For 
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ordered systems on the other hand, polarons (20; 21; 22; 23; 24) result in ad- 
ditional peaks in the Hubbard bands. As was noted by Sangiovanni et al. (25) 
these polaron peaks are present in dynamical mean field theory (26) but not in 
the static LDA+U mean field theory. Another drawback is that SIC-LDA and 
LDA+U almost automatically yield localised electrons with split (Hubbard) 
bands, even if this is not correct. 

Missing in both LDA and LDA+U is also the strongly correlated metallic 
phase found in many transition metal oxides and heavy Fermion systems either 
at intermediate values of C/ or at a non-integer number of d or / electrons per 
site so that the Mott insulator is doped and becomes metallic. This strongly 
correlated metal already has Hubbard bands as the Mott insulator, but the 
low-energy behaviour is dominated by quasiparticle physics, see the central 
panel of Fig. 2. We can think of these quasiparticles as dressed electrons which 
move independently from each other, although with a larger effective mass m 
since the quasiparticle consists of the initial electron plus the electron-electron 
interaction with its environment. The mass enhancement m/nie becomes more 
pronounced with increasing U, and the weight of the quasiparticle peak Z = 
{m/mf,)~^ is reduced correspondingly until it vanishes completely at the Mott- 
Hubbard transition where the system becomes insulating. 

A modern, non-perturbative technique to describe these many-body effects 
is the dynamical mean field theory (DMFT) developed primarily by Metzner 
and Vollhardt (27), Miiller-Hartmann (28), Georges and Kotliar (29), and 
Jarrell (30); for more details see Section 3. One of the strong points of DMFT 
is that it can describe such strongly correlated metals and Mott insulators in 
a single framework. 

Recently, physicists from the many-body and DFT community have devel- 
oped a new approach to correlated materials which merges two of the most 
successful approaches of the two communities: LDA-I-DMFT. The initial work 
by Anisimov et al. (31) and Lichtenstcin and Katsnclson (32) was followed by 
a rapid development which is the subject of this Review. Shorter LDA-I-DMFT 
summaries have already been published in the form of conference proceedings, 
lecture notes and a Psi-k Highlight, sec Rcfs. (33; 34; 35; 36; 37; 38; 39); a full 
review was submitted by Kotliar et al. (40) during the completion of this work. 
Instead of LDA, one can also use the so-called GW approximation of Hedin 
(41), resulting in the related GW-I-DMFT approach formulated by Biermann 
et al. (42) and presented in more detail in Ref. (43) and Ref. (44). 

The LDA+DMFT approach does not only contain the correct quasiparticle 
physics and energy of the strongly correlated metallic phase, but also repro- 
duces the correct results in the limit of small and large C/, see Fig. 2. Hence, 
LDA+DMFT correctly accounts for the main contributions of electronic cor- 
relations (there are corrections due to the k-dependence of the self energy as 
discussed in Section 3.2). In contrast, LDA yields an uncorrelated metal even 
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for strongly-correlated metals or Mott insulators. LDA+U on the other hand 
typically predicts an insulator for the ab-initio calculated U values of 3d tran- 
sition metal oxides, even for materials which are metallic. Considering this, 
it might not be so astonishing that LDA-I-DMFT turned out to be a break- 
through for the calculation of materials with strong electronic correlations. 

In the following, we will review these recent developments at the former bor- 
derline between DFT and many-body theory. Thereby we take a Hamiltonian 
point of view which allows us to understand physically what LDA-I-DMFT 
describes, and what it does not. Complementarily, the review of Kotliar et al. 
(40) focuses on the spectral density functional theory for which LDA-I-DMFT 
is the standard approximation just as LDA is for DFT. 

We start, in Section 2, with a brief recapitulation of conventional DFT/LDA 
band structure calculation. For more LDA and DFT details, we advise the 
reader to consult the excellent review by Jones and Gunnarsson (3) or the 
books Refs. (3; 4; 5; 6). In the present Review, we will not elaborate on 
how to construct an optimal set of wave functions within LDA and how to 
choose maximally localised orbitals for the further many-body DMFT calcu- 
lation. To this end, the linearised muffin tin orbitals (LMTO) of Andersen 
(45), for example, in the third generation (46; 47), or the Wannier (48) func- 
tion projection by Marzari and Vanderbilt (49) can be employed; also see Refs. 
(50; 51; 52; 53; 54; 55). We refer the reader to these original and review pub- 
lications for the full specification of the orbital basis sets whose construction 
is a science on its own. 

In Section 3, we give a compendious introduction to DMFT, without dis- 
cussing all derivations of the DMFT equations and the physics of the Mott 
transition since this has been reviewed by Georges et al. (56). The central part 
of this review is Section 4 in which the recently proposed methods employing 
DMFT for realistic material calculations are introduced. In particular, we will 
discuss the LDA-I-DMFT method in Section 4.1, a combination of the Hartree 
approximation and DMFT in Section 4.2, and the GW-j-DMFT approach in 
Section 4.3. For more details on the GW approximation itself, we recommend 
the review by Aryasetiawan and Gunnarsson (57). Section 5 is devoted to sev- 
eral methods, numerically exact and approximate ones, for solving the DMFT 
equations. The respective advantages and disadvantages are discussed. 

In Section 6, we review the successful electronic structure calculations with 
DMFT accomplished for /-electron systems, transition metals and their oxides. 
Finally in Section 7 we give a summary and outlook. 



February 2, 2008 



18:28 Advances in Physics ESCwDMFT'16 



8 K. Held 

2 Conventional electronic structure calculations 

2.1 Ab initio Hamiltonian 

Let us start by rewriting the ah initio Hamiltonian [1] in second, instead of 
first, quantised form: 



47reo |r-R; 

syi„„(r) 
2 



+ \Y. jd''rd\'i>+{v,a)i>+{v',a') *(r', a')^(r, C7),(2) 



=^^ee(r-r') 



where now ^"'"(r, o") and ^(r, o") arc field operators that create and annihi- 
late an electron at position r with spin a; Vion(r) and Vcd^) denote the lattice 
potential and the electron-electron interaction, respectively; for the other sym- 
bols see Eq. [1]. 

Not included in Eq. [2] are the kinetic energy of the lattice ions (phonons) 
and relativistic corrections, in particular, the spin-orbit coupling. Both are 
known to be important in some materials. For example, the electron-phonon 
coupling plays an important role in Jahn- Teller distorted manganites according 
to Millis et al. (58). The spin-orbit coupling on the other hand is important 
for /-electron systems. 



2.2 Density functional theory (DFT) 

The basic idea of density functional theory (DFT) is to work with a simple 

quantity, i.e. the electron density p(r), instead of trying to solve the ah ini- 
tio Hamiltonian [1,2] through complicated many-body wave functions. This is 
possible, at least for the ground state energy and its derivatives, thanks to the 
Hohenberg-Kohn (2) theorem, which states that the ground state energy is a 
functional of the electron density £^[p(r)] which is minimised at the ground 
state density p(r) = po{v). Following Levy (59), this theorem can be easily 
proven and the functional even be constructed by a two-step minimisation 
process, see Fig. 3. In a first step, we select all many-body wave functions 
V3(riO"i, ... rATCJAr) for a fixed number of electrons N which yield a certain elec- 
tron density p(r), and minimise the energy expectation value for these wave 
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all wave functions 



wave functions with p 




wave functions with p. 



Figure 3. Construction of the DFT functional according to Levy (59). The DFT energy functional 
for a given electron density pi (r) is the minimal expectation value of all wave functions yielding 

pi{r). 



functions 



E[p] 



mm 



N 



(3) 



This is exactly the DFT functional for the ground state energy: A second min- 
imisation yields the ground state energy Eq = miiip E[p] at the ground state 
density po, since for the ground state density E[pq] includes the energy expec- 
tation value w.r.t the ground state wave function, see Fig. 3 for an illustration. 

While this construction proves the Hohenberg-Kohn theorem, we did actu- 
ally not gain anything: For obtaining E[p] we have to calculate the expectation 
value {(p\H\ip) for complicated many-body wave functions Only the ionic 
potential Eion[p] = / d^r Mon(r) p(r) and the Hartree term i^Hartree [/o] = 
^ J d?r' d?r 14e(r — r') can be expressed easily through the electron 

density. Denoting by -E'kin[/'] the kinetic energy functional, we can hence write 



E[p] = -Ekinl/O] + Eion[p] + -EHartreeM + Ey^c[p\, 



(4) 



where all the difficulty is now hidden in the exchange and correlation term 
E^c[p\- This term is unknown. An important aspect of DFT is however that 
the functional E[p\ — E\o^[p] does not depend on the material investigated; 
for a prove see e.g. Ref. (3). Hence, if we knew the DFT functional for one 
material, we could calculate all materials by simply adding E\oa[p\- 

For calculating the ground state energy and density, we have to minimise 
5{E[p\ - p{J d^rp{r) - N)}/dp{r) = where the La grange parameter p fixes 
the number of electrons to N. To avoid the difficulty of expressing the kinetic 
energy E'kinl/o] through p{r), Kohn and Sham (60) introduced an (at this point 
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auxiliary) set of one-particle wave functions cpi yielding the density 

N 



p{r) = Y,\^,{r)\\ (5) 



and minimised w.r.t. the (/?.j's instead of p(r). That is, wc minimise — 
EilJ (i^r|(/?j(r)p] — 1}/Sipi{r) = where the Lagrange parameters guarantee 
the normalisation of the ipi's. This minimisation leads us to the Kohn-Sham 
(60) equations 



^^-A + M„.(r)+/d3/Fee(r-r')p(r')+ ^^^^^ 



ipi{r) = £i (^j(r). (6) 



These Kohn-Sham equations are Schrodinger equations, describing single elec- 
trons moving in a time-averaged potential 

Kff(r) = Mon(r) + I d\' Fee(r-r')p(r') + (7) 

of all electrons. The Kohn-Sham equations and the electron density have to 
be calculated self-consistcntly, see the flow diagram Fig. 4. 

Let us note that the one-particle Kohn-Sham equations serve, in principle, 
only the purpose of minimising the DFT energy, and have no physical meaning. 
If wc knew the exact E'xc, which is non-local in p(r), wc would obtain the 
exact ground state energy and density. However, in practise, one has to make 
approximations to -Exc such as the local density approximation discussed in the 
next Section. We can then think of these approximations as describing single 
electrons moving in an approximated potential 6Eyj^\p\lbp{x)^ as illustrated in 
the Introduction, Fig. 1. 

Let us also note that the kinetic energy in Eq. [7], -EkiniPmin] = 
— ^^^((/?i|/i^A/(2me)|</'i), is that of independent (uncorrelated) electrons. 
The true kinetic energy functional for the many-body problem is different. 
We hence have to add the difference between the true kinetic energy func- 
tional for the many-body problem and the above uncorrelated kinetic energy 
to E'xc, so that all many-body difficulties are buried in E'xc- 



2.3 Local density approximation (LDA) 

Since we do not know E^c exactly, wc have to make approximations, and 
the most widely employed approximation is the local density approximation 
(LDA). In LDA, the complicated non-local functional -Excfp] is replaced by a 
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Choose initial electronic density p(r) 




Solve Kohn-Sham equations [6]: 

["f ^ + ^on(r)+ J dh%,ir rOp(r') + ^^ffy^ir) = sm{r) 


Calculate new electronic density via Eq. [5] : 

AT 

Pnew(r) = ^\(pi{r)f. 

i 


Iterate with p{r) = /9new(r) until convergence, i.e. ||p(r) — pnew(r)|| <e. 



Figure 4. Flow diagram of the DFT/LDA calculations. 

local (LDA) exchange energy density^ which is a function of the local density 
only: 

Jd'rE^^\p{r)). (8) 

In practise, E^^^{p{r)) is calculated from the perturbative solution (61; 62) 
or the numerical simulation (63) of the jellium model which is defined by 
Mon(r) = const. Due to translational symmetry, the jellium model has a con- 
stant electron density p{r) = pQ. Hence, with the correct jellium E^^, we 
could calculate the energy of any material with a constant electron density 
exactly. However, for real materials p(r) is varying, less so for s and p valence 
electrons but strongly for d and / electrons. While LDA turned out to be 
unexpectedly successful for many materials, it fails for materials with strong 
electronic correlations between the d or / electrons. In this paper, we will 
review recent advances to improve on LDA for such correlated materials by 
using DMFT. 

To understand these advances physically, let us now turn to the LDA calcu- 
lation of band structures, one of the main applications of LDA. If we calculate 
a LDA band structures we will leave the firm ground of DFT, which strictly 



^Notc that in the literature, e.g. in (3), often the exchange energy density per particle e^^'^(p(r)) = 
ELDA(p(r))/p(r) is used. 
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speaking only allows for calculating the ground state energy and its derivatives. 
Instead, for the LDA bandstructure, we interprete the Lagrange parameters 
as the physical (one-particle) energies of the system. This corresponds exactly 
to the picture Fig. 1 of Section 1. Physically, the one-particle LDA bandstruc- 
ture corresponds to approximating the ah initio many-body Hamiltonian by 



HLDA=Eld^r ^+(r,a) 



^2 



2m, 



■A + Vion(r)- 



+ 



dpir) 



*(r,a).(9) 



In practise, one solves the three-dimensional one-electron Kohn-Sham equa- 
tions given by Hamiltonian [9] through expanding the wave functions in sophis- 
ticatedly chosen basis sets. These basis sets use either the atom and so called 
muffin tin potentials, or free electrons in, e.g. ultrasoft (64), pseudopotentials 
as a starting point; for a detailed discussion sec (4). An important aspect for 
merging LDA with DMFT is that wc need to identify the correlated d or f 
orbitals. For the linearised muffin-tin orbital (LMTO) basis of Andersen (45), 
and its Nth order extension (NMTO) (46; 47), these orbitals are inherently 
defined. For other basis sets, e.g. plane waves, one has first to project onto 
Wannier orbitals (49) before doing the DMFT calculation. 

Let us now transform our r-basis to such a basis with localised d or / orbitals: 



*+(r,c7) = ^c[M,(r). 



(10) 



Here, I denotes the orbital index (which, for problems with more than one 
atom in the unit cell, will also index orbitals on different sites in this unit 
cell.) Transforming the Hamiltonian [9] to this orbital basis, we obtain 



-^LDA — ^ tiljm (^il^^jm^ 
il jm a 



(11) 



where tujm denotes the tight-binding hopping matrix elements 



til 



jm 



2mp 



+ T4on(r) + 



dV'p(rOKe(r-rO + 



dE^^^{p{v)) 

dp{v) 



<^,-m(r). (12) 



According to the Bloch theorem, we can diagonalise Hamiltonian [11] by 
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Fourier transforming to Ck; = J2i (kl^~^^^ {L denotes the number of lattice 
sites Ri): 

^LDA = ^^^-^(k) C^lc^^. (13) 

Here, 

e^r(k) = ^E*^Ome^'^(^^-^^-) (14) 

is still a matrix in our orbital indices which can be fully diagonaliced by a sim- 
ple orbital rotation. This orbital rotation yields the Kohn-Sham eigenvalues, 
i.e. the LDA bandstructure. Let us emphasise that this rotation is different 
one for every k point. 

As already mentioned, LDA is highly successful for calculating both, static 
physical properties and bandstructures. But LDA is not reliable when applied 
to correlated materials and can even be completely wrong because it treats 
electronic correlations only very rudimentarily: The exchange correlation func- 
tional has been approximated by a functional which depends only on the local 
density. It has been calculated from the jcUium problem, a weakly correlated 
model with extended orbitals. Basically, LDA is a one-electron approach where 
the effect of the other electrons is through a static (time independent) mean 
field. In the following, we will discuss how to treat electronic correlations in a 
better way with dynamical mean field theory (DMFT). 



3 Dynamical mean field theory (DMFT) 

The conventional approach to strongly correlated electron systems is perturba- 
tion theory, either starting from the limit of zero Coulomb interaction {U = 0; 
weak coupling) or from the opposite limit of zero kinetic energy {W = 0; 
strong coupling). Such a treatment is, however, not appropriate for the im- 
portant 'in-between' regime where we have a strongly correlated metal and 
possibly a Mott-Hubbard transition between this correlated metal and a Mott 
insulator (sec Fig. 2). 

Metzner and Vollhardt (27) introduced a new limit to correlated electron 
systems, the limit of infinite dimensions d — ^ oo or equivalently an infinite num- 
ber of neighbouring lattice sites. In this limit, the competition between kinetic 
energy W and Coulomb interaction U is maintained, albeit resulting in a sim- 
plified, momentum- independent self energy: Sk(t(;) ''-I^ S(ti;). The work by 
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Metzner and Vollhardt (27) initiated a rapid development: Miiller-Hartmann 
(28; 65; 66) showed only the local Coulomb interaction yields dynamic {uj- 
dependent) correlations whereas the non-local density-density interactions are 
reduced to the Hartrcc contribution which is w-indcpcndcnt. Some models like 
the Falikov-Kimball model (67) and the Kondo lattice model for classical spins 
could even be solved exactly by Brandt and Mielsch (68; 69; 70; 71; 72; 73) 
and Furukawa (74; 75), respectively. The next important step was to put the 
initial ideas on the footings of a mean field theory, coined dynamical mean 
field theory (DMFT) because of the frequency dependence of E(a;). For the 
Falicov-Kimball model, Brandt and Mielsch (68; 69; 70) already succeeded in 
this respect, and Janis and Vollhardt (76; 77) generalised the coherent poten- 
tial approximation (CPA) to include the dynamics of correlated electron mod- 
els in the d ^ oo limit. But the breakthrough represents the work by Georges 
and Kotliar (29) (also see Ref. (30) and Ref. (78; 79)) who showed that a 
many-body model like the Hubbard model is mapped onto the self-consistent 
solution of an auxiliary Anderson impurity model for d oo (for the physics 
of this Anderson impurity model see e.g. (80; 81)). This mapping was crucial 
for the following development since physicists were thenceforth able to employ 
well-known solvers for the Anderson impurity model to deal now with lattice 
many-body models in the DMFT framework: Approximate solvers like iterated 
perturbation theory (29), self-consistent perturbation theory (66; 82), and the 
non-crossing approximation (83; 84; 85; 86); and numerically exact solvers like 
quantum Monte Carlo simulations (30; 87; 88), exact diagonalisation (89) and 
the numerical renormalisation group (90; 91). 

The DMFT results turned out to be a big step forward for our understand- 
ing of, for example, the Hubbard model and the periodic Anderson model, 
and improved our insight, particularly into the Mott-Hubbard transition; see 
Georges et al. (56) for a review and Kotliar and Vollhardt (92) for a first 
reading. 

Since DMFT is central to this review, which addresses readers from different 
communities, we feel a self-contained derivation of DMFT is helpful. In Section 
3.1, we present such a derivation based on the local nature of the weak-coupling 
perturbation theory, as first noticed by Metzner and Vollhardt (27) and Miiller- 
Hartmann (28). Alternative derivations of the DMFT equations have used 
the cavity method (56), the expansion around the atomic limit (93) and a 
generalisation of the coherent potential approximation (76: 77). The reader 
who is not interested in this derivation might directly turn to the DMFT self- 
consistency scheme summarised on p. 21 ff. Shortcomings and extensions of 
DMFT by non-local correlations are discussed in Section 3.2. 
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3.1 Derivation of the DMFT equations 

Let us in the following consider a generalised many-body Hamiltonian 

^= ^iljm cjcj^+ ^ Uimno C-ii^ C-jJ clj^clg, (15) 

il jm a i Imno era' 

where c^^ (c^i) creates (annihilates) an electron with spin a and orbital index 
I at lattice site i; tnjm is a hopping amplitude between lattice sites i and j 
and orbitals I and m [of the same form as in Eq. [11]); finally, Uimno denotes 
a general local Coulomb interaction. 

Proper scaling in the limit d ^ oo- Let us consider an extension of the lattice 
at hand so that we have a large number Zj^^j of equivalent (neighbouring) 
sites, with the same 'distance' j|| to site i. Later, we will consider the limit 
GO. The reader might envisage, for example, the generalisation of 
the cubic lattice to a d-dimensional hypercubic lattice. This hypercubic lattice 
has Zi = 2d (6 for d = 3) nearest neighbours and, in general, ~ di*"^!! 

equivalent sites at distance ||?' — jI (to leading order in 1/d). For the cubic 
lattice, II ... II is the so-called Manhattan norm. The advantage of the notation 
is that all classes of equivalent sites and different lattice topologies can 
be treated on the same footings. The reader will however realise that these 
Zjj_jjS are in general not independent. 

How do the two terms of Hamiltonian [15] scale with increasing Zj^^jj in 
the limit Zjj^jj oo7 Obviously, the purely local interaction Ui^nno and the 
associated potential energy per site scales like 

( Yl UlrnnoCl^C^^C^nCfo I const., (16) 

\lmnoaa' I 

i.e. stays at a finite constant, neither going to zero nor infinity. Here and in 
the following, 

Tr de~^^ 

denotes the thermal expectation value at inverse temperature /3 = l/T; Tr is 
the trace. 

For the first term of Hamiltonian [15], the kinetic energy, we have Zjj_^j 
equivalent j terms in the sum. Hence, the kinetic energy per site i diverges 
unless we rescale tujm in the limit oo. To avoid this divergence, let 
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US consider the following ansatz for the scaling of these hopping amplitudes: 

f* 

tiljm - ,„ UoJ 

V^IIHII 

with t*i^^ staying constant upon increasing -S^ii-ji- Since the non-interacting 
{U = 0) Green function^ G^ij^{u) is directly connected to tujm it scales 
exactly in the same way^: 

G^i,M = [{col-t)-'] ~ -y^. (19) 

Here, bold symbols with double underlines denote matrices w.r.t. orbital and 
site indices {1 is the unit and t the hopping matrix) . Let us assume for a mo- 
ment that the very same scaling also holds for the interacting Green function 

GujM ~ -T^- (20) 

Later, we will see that indeed all Feynman diagrams for G{lo) scale like Eq. 
[20] or fall off even faster. 

With the scaling Eq. [20], the kinetic energy per site scales properly, i.e. 
stays finite: 

Yl ^^^i^ i^a^jm) ^'"^"^ const. (21) 

IjmcT 



To see this, note that there are terms in the sum for every class of 

]_ 



equivalent sites j. This factor is however canceled by a factor —7= — for Uijm 



^Depending on the problem at hand, we will work with Green functions for real timc/frcqucncics or 
imaginary time/Matsubara frequencies [uii, = {21/ + l)7r//3]; see Ref. (94) for an introduction. Both 
representations are connected via the analytical continuation in the complex plane. Let us define for 
completeness the real time Green function 

Giljmit - t') = -iQit - t'){Tc.fi)c]^{f)) 

and the imaginary time (thermal) Green function 

Giijmir - r') = -{Tc,i{-iT)cl^{-iT')). 

Here, T is the Wick time-ordering operator w.r.t. t or r and Cj;(t) = exp(i/ft)cj; cxp(—iHt) in both 
cases. The corresponding frequency-dependent Green functions follow by Fourier transformation. 
^Expressing the matrix elements of the inversion in Eq. [19] in terms of minors one directly sees that 
the off-diagonal elements are a factor tij j,„/(w — tjiim) smaller than the diagonal elements. 
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Figure 5. A second-order diagram for the Green function (or the self energy if the legs are 
amputated). A simple power counting of the three Green functions connecting i and j shows that 
the contribution of this diagram vanishes in the d — > oo limit. Only the local contribution with 

i = j survives. 

and another factor , ^ for {(fj (fj^) which is directly connected to the 
Green function 

Gujm{r = 0+) = -5ij5im + {c-!c'^J. (22) 

Only the scaling of Eq. [18] provides a non-trivial limit in which the interplay 
between kinetic and Coulomb interaction energy remains intact. 

Locality of diagrams- Now, let us consider the Feynman diagram for the 
Green function depicted in Fig. 5. The contribution of this diagram, but gen- 
erally every diagram where two sites i ^ j are connected by three independent 
lines (G^j's)^ scales like l/^yZj~~j^ for individual i's and j's. This is a factor 
1/Zji_jj smaller than the direct contribution Eq. [19], i.e., a Green function line 
without interaction. Hence, such non-local diagrams for the Green function (or 
self energy) become irrelevant in the limit oo. 

There are also Feynman diagrams as in Fig. 6, where i ^ j are connected 
by only two lines. Taking into account a factor for the sum over differ- 

ent j, we see that the overall contribution of this diagram scale like Zjj^jj'. 
However, if we write the Feynman diagrams in terms of the interacting (full) 
Green function Gij instead of the non-interacting G^^ this diagram is already 
contained in the local diagram on the right hand side of Fig. 6. This is because 
Gij includes self-energy inclusions as on the left hand side of Fig. 6. 

Generally, the relation between Gij (double lines) and G^^- (single lines) is 
given by the Dyson equation which reads 

G,,(r) = G^,ir) + Y.f^r'j dr" Gl,{T')^,r{r" - r')G,.,(T - r") (23) 






^Here and in the following, bold symbols denote matrices in the orbital index. 
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Figure 6. Left hand side: A diagram for the Green function or self energy which contributes in the 
d — > oo limit since i and j are connected only by two Green function lines. However, this diagram is 

contained in the local diagram on the right hand side in terms of the full (interacting) Green 
function (double line), since it includes interaction or self energy parts as the one on the left side. 
For d — > oo, all {skeleton) diagrams in terms of the interacting Green function are purely local, i.e. 

only involve one site i. 



or graphically 




Fourier transformed, Eq. [23] simplifies to 



Gk(iw,)-^ = Gliiuj,)-^ - ^^{iuJu). (24) 

If we now write our diagrams in terms of Gjj, we have to restrict ourselves 
to a subset of all Feynman diagrams, the so-called skeleton diagrams which 
do not contain any part connected to the rest of the diagram by only two Gjj 
lines because such diagrams are generated if Gij is expanded in terms of G^^ . 
If these were included in the skeleton diagrams they would we counted twice. 

How do the skeleton diagrams scale with If there are two sites i ^ j, 

both sites have to be connected by three Green function lines (directly or 
indirectly via additional sites), which follows directly from the definition of 
the skeleton diagrams. Therefore, the scaling for individual sites j ^ i goes at 
least like Xj ^JZ^y^^ as in Fig. 5. These contributions become irrelevant for 
Z\i-j\ oo. Consequently, all skeleton diagrams are purely local. An electron 
may still leave site i, interact on other sites, and return to site i retardedly. But 
these processes are entirely described by the local interacting Green function 
Gii, containing self energy inclusions as in Fig. 6. 

As a result, the self energy itself is purely local: 



(25) 
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Figure 7. DMFT maps the lattice many-body problem with interactions U on every site (left side) 
onto a single site problem where the interaction has been replaced by the self energy S except for a 
single site (right side). All irreducible Feynman diagrams are hence local: electrons, leaving and 
returning to this single interacting site, are dressed by S (double line in Fig. 6). The DMFT 
mapping becomes exact for d —* oo and is an approximation in finite dimensions. 



or Fourier-transformed k- independent: 5]k(a;) = 1]{uj). From this and the 
Dyson equation [23], it follows that Gujmi^) ~ ^/ ^ posteriori con- 
firming our assumption Eq. [20]. 

To conclude this scaling analysis, we can say that DMFT represents the 
local contribution of all skeleton Feynman diagrams, thereby replacing the 
lattice many-body problem by a local (single site) problem, see Fig. 7. This 
simplification becomes exact for d — > oo and is an approximation if applied to 
a finite dimensional problem. 

Mapping onto the Anderson impurity model- Diagrammatically, DMFT cor- 
responds to the local contribution of all topologically distinct Feynman dia- 
grams. Exactly the same diagrams can be obtained via an Anderson impurity 
model if its on-site interaction has the same form as the original Hamiltonian 
[15]: 



Imno aa' 



Here, (^kz) creation and annihilation operators for non-interacting con- 
duction electrons at wave vector k which have a dispersion ei(k) and hybridise 
with the localised interacting electrons (fm via V;m(k). 

Let us now switch to the language of functional integrals in terms of Grass- 
mann variables 'i/' and "0*) see e.g. Ref. (95). By a simple Gaussian integration 
we can get rid of the conduction electrons since these only enter quadratically 
in Hamiltonian [26] and, hence, in the functional integral. Integrating out these 
conduction electrons, we arrive at an effective problem for the interacting elec- 
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trons Cm^ whose Green function reads 

Here, v denotes imaginary Matsubara frequencies Ui, = ■k{2u + 1)//?; 

Z = J D[V']X)[V'*]e'^['^'^*'(^'')"'l (28) 
is the partition function and the single-site action A has the form 



ualm 

/3 



- E ^irnnoj dr ^f* {r)^n {r)^^* {rWo (r). (29) 

Imno aa' q 

In Eq. [29], the interaction part of A is in terms of r, the Fourier transform 
of the Matsubara frequencies Ui^. The non- interacting Green function of the 
effective problem is given by 

[GZi^^^r' = + Un.rn + 1^-^! .^"ff^^^^tV ^^^^ 

^ lUv + fi — ei[k) 

with the same local term timin as the original many body Hamiltonian [15]. 

The topology of the irreducible diagrams of this effective Anderson impurity 
model is exactly the same as the DMFT single site problem: simply the local 
contribution of all Fcynman diagrams. If also the lines of the diagrams arc the 
same we will hence obtain the same self energy. For the local irreducible dia- 
grams, this is achieved if the interacting Green function G(a;) of the Anderson 
impurity model also equals the local DMFT Green function Gii{uj). 

Which Anderson impurity model has now this G{uj) as its interacting Green 
function? -"^ Given its interacting Green function G(w) and self energy S(w), 
the impurity model's non-interacting Green function has to fulfil 

[G\u)]-' = [G{ur' + ^{u). (31) 



^ Since Green function G(uj) and self energy ^{ui) are the same for DMFT and the Anderson impurity 
model we employ the same symbols. 
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This is the analogy to the Dyson equation [23], but now for the auxiUary An- 
derson impurity model. Note, that if we expand the G{io) of the irreducible 
diagrams in terms of G^{i^), we will get reducible inclusions and all local Feyn- 
man diagrams. While the reducible diagrams for the Anderson impurity model 
are different from the (non-local) reducible DMFT diagrams, the irreducible 
diagrams are exactly the same. Hence, we can calculate the interacting DMFT 
Green function or self energy by solving an Anderson impurity model. 

Self consistency scheme- Building upon these considerations, we can now 
formulate the self-consistent DMFT scheme, which is summarised in the 
DMFT flow diagram Fig. 8. We start with a trial self energy S(a;), which 
might just be zero. Then, we calculate the local Green function G(a;) by k- 
integration of the Dyson equation [24] : 



Here, we have made use of G°(k)-i = Lol + fil - e(k) where e(k) is the 
Fourier transform of tujm like in Eq. [14], 1 the orbital unit matrix, and Vbz 
the volume of the Brillouin zone (BZ). 

Given G(a-') and XI (oj), we define the effective non-interacting Green function 
of the Anderson model G^{i^) via Eq. [31]. Now comes the difficult part: We 
have to solve the Anderson impurity model and calculate its Green function 
G(a;). Different methods to this end are presented in Section 5. Using Eq. 
[31] a second time, we obtain a new self energy. (Sometimes, in particular in 
perturbation theory, it is more convenient to calculate the self energy of the 
Anderson impurity model directly, like in our diagrammatic presentation in 
Section 3.1.) With this new self energy 5](a;), the algorithm is iterated until a 
convergence criterion e in some norm ||...||, e.g. the maximum or sum norm, is 
met. 



3.2 Extensions of DMFT 

DMFT accounts for a major, local part of electronic correlations. This part 
is responsible for the strong quasiparticle renormalisation and the Mott- 
Hubbard transition in materials with strong electronic correlations such as 
transition metal oxides and heavy Fcrmion systems, see the discussion in 
Section 1. But also many other physical phenomena have been described 
by DMFT: the crossover from Slater (96) to Heisenberg (97) antiferromag- 
netism (30) and associated metamagnetism (98; 99), the infiucnce of the lattice 
(100; 101) and Hund's exchange on ferromagnetism (102), scattering at impuri- 
ties (103; 104; 105; 106) or spins (74; 75), the effect of phonons (107; 108; 109), 




(32) 



BZ 
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Choose an initial self energy S. 




Calculate G from Xl via the k-integrated Dyson equation [32] : 

G(a;) = / d^k[ujl + iil- e(k) - ^{uj)]~^ . 
Vbz J 

BZ 


e?o = (G-1 + (Eq. [31]). 


Calculate G from via the DMFT single-site problem Eq. [27]: 


5]new = - (Eq. [31]). 


Iterate with 5] = Xlnew until convergence, i.e. — 5]new|| < 



Figure 8. Flow diagram of the DMFT self-consistency cycle, also see Ref. (35). 

particularly on the 'colossal magnetoresistance' in manganites (110; 111; 112) 
and superconductivity (113; 114). This list is far from complete, with only a 
bare minimum of references. Nonetheless, it shows that DMFT yields a wide 
variety of physical phenomena, posing the question: What does DMFT not 
describe? 

What DMFT neglects completely are non-local correlations. Such non-local 

correlations are responsible for a variety of physical phenomena, ranging from 
valence bonds, pscudo gaps and (possibly) d-wave superconductivity to (para- 
)magnons, quantum critical behaviour, and generally the critical behaviour in 
the vicinity of phase transitions. All these effects are not included in DMFT. 
Some of these effects stem from rather short range correlations, e.g. the va- 
lence bond formation of pairs of spins into singlets; others require long-range 
correlations, e.g quantum criticality. 

There have been various attempts to include non-local correlations be- 
yond DMFT. Most successfully and widely employed are cluster extensions 
of DMFT. Since these have been reviewed recently by Maier et al. (115), we 
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Figure 9. In cluster DMFT the lattice many-body problem (left side) is mapped onto a DMFT 
cluster (right hand side) which consists of several (here two) interacting sites. Nonlocal correlations 
between the cluster sites are taken into account. The DMFT field couples to the border of the 

cluster. 



will only briefly discuss the basic ideas and different cluster schemes, without 
giving referenced to all applications to model systems. 

The idea of these cluster extensions is to treat, instead of the single DMFT 
site, a whole cluster of sites, see Fig. 9. A natural choice for the cluster is 
a super cell in real space, or the LDA unit cell if it contains more than one 
site with interacting d or / orbitals. Then the DMFT basis block is not a 
single site but this super cell. Non-local correlations within the super cell are 
taken into account, whereas such correlations between different super cells 
are neglected. Beyond the local DMFT self energy, there are now also off- 
diagonal elements 'Sij{Lo) for two sites i ^ j within the cluster. This path is 
followed in the so-called cluster DMFT approaches (56; 116; 117; 118). Cluster 
DMFT turned out to yield surprisingly precise results for the one-dimensional 
Hubbard model (119). With the number of cluster-sites, the numerical effort 
grows considerably, in particular if realistic calculations require the inclusion 
of several orbitals. Such realistic LDA+cluster DMFT calculation have been 
employed so-far by Poteryaev et al. (120) for studying Ti203, by Mazurenko 
et al. (121) for NaV205, and by Biermann et al. (122), who found VO2 to be 
insulating because spins form non-local singlets so that a Peierls gap opens. 
In both cases the LDA unit cell provided for a natural choice of a two-site 
cluster. 

Another scheme, named dynamical cluster approximation (DC A), has been 
proposed by Hettler et al. (123). While cluster DMFT can be best understood 
in real space, the idea of DCA is to patch the self energy in k-space. To this 
end, the first Brillouin zone is divided into N^. patches around k- vectors K. 
In contrast to the constant DMFT self energy, the DCA self energy Sk is 
only constant within the patch, but varying from patch to patch. This way, 
k-dependencies of the self energy are taken into account. The basic differ- 
ences between DCA and cluster DMFT in real space are: The DCA cluster 
has periodic boundary conditions, whereas that of cluster DMFT has open 
boundary conditions. Moreover, the dynamical mean field couples to every 
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site of the cluster in DCA and to the boundary sites in cluster DMFT. In the 
limit of infinite cluster size (A'^c oo), both approaches become exact. While 
cluster DMFT might be more suitable for dealing with particular short-range 
correlations in real space, such as the formation of a spin singlet on neighbour- 
ing sites, DCA more naturally preserves translational symmetry. Let us also 
mention three further alternatives: the cluster perturbation theory by Gros 
and Valenti (124), the variational cluster perturbation theory of Potthoff et 
al. (118), and the (very general) self energy functional approach of Potthoff 
(125; 126). In the context of classical spin models cluster extensions such as 
the embedded cluster method (127) have been used for a long time, leading to 
early proposals (56) to employ similar approaches for DMFT. 

In the DMFT context, cluster extensions have been applied extensively by- 
now for studying superconductivity. Indeed, cZ-wave superconductivity in the 
two-dimensional Hubbard model was found for a 2 x 2 cluster (128; 117; 129), 
and a pseudogap phase at elevated temperatures (130). Whether the d-wave 
superconductivity persists for larger cluster sizes and in the Nc oo limit is 
however less obvious (115; 131; 132). Let us also note that Jarrell and coworkers 
(133; 134) employed a hybrid cluster solve for treating larger cluster sizes, 
where the the fluctuation exchange approximation (FLEX) by Bickers et al. 
(135) is employed for longer-range correlations. 

A different route to extend DMFT by non-local correlations has been pro- 
posed by Schiller and Ingersent (136; 137). The authors extended DMFT dia- 
grammatically through all diagrams to order 1/d. This leads to a theory with a 
single-site and a two-site cluster whose Green functions have to be subtracted. 
Necessarily, this 1/d approach is hence restricted to nearest neighbour corre- 
lations. 

All extensions of DMFT mentioned so-far, are restricted to relatively short- 
range correlations. Often, however, long-range correlations are of vital im- 
portance, e.g. for magnons, screening and quantum criticality. If we want to 
account for such long-range correlations, other extensions of DMFT are nec- 
essary. Very recently, Toschi et al. (138) introduced a way to do so (also see 
Ref. (139; 140)): the dynamical vertex approximation (DFA). The basic idea 
is to extend DMFT diagrammatically: Toschi et al. (138) take the local, fully 
irreducible two-particle vertex as a starting point and construct from this ver- 
tex all possible (local and non-local) self energy diagrams. A restriction of 
this approach to the particle-hole channel yields (para-)magnons, a restriction 
to the particle-particle channels yields the cooperon diagrams. Hence, ladder 
diagrams describing long-range correlations for weakly correlated systems are 
recovered, but with the vertex instead of the bare Coulomb interaction so that 
strong correlations are accounted for. Independently, Kusunose (139) extended 
DMFT by the transverse particle-hole channel using iterated perturbation the- 
ory, and Slezak et al. (140) proposed to use similar ladder diagrams for solving 
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larger DCA clusters. Previous attempts supplement the DMFT self energy by 
spin-fluctuations from the spin-fermion model (141; 142; 143; 144) or from the 
self-consistent renormalisation theory (145). Bolech et al. (146) explored the 
possibility of a renormalisation group extension of cluster DMFT and Rubtsov 
et al. extended DMFT by including non-local correlations through the dual 
Fermion approach (147; 148). 



4 Merging conventional bandstructure approaches with DMFT 

The standard 20th century approach to calculate materials realistically, i.e. 

density functional theory (DFT) in its local density approximation (LDA), 
was introduced in Section 2. As pointed out in Section 1, this conventional 
approximation does not work for materials with d or f electrons because elec- 
tronic correlations arc too strong. Of prime importance for these correlations 
are the local Coulomb interactions between the d ox f electrons. This is be- 
cause the local Coulomb interaction is the largest one and also because other 
density-density interactions are to leading order in l/Z^j^^ (one over the num- 
ber of sites with the same distance, sec Section 3) only given by the Hartree 
term as was proved by Miiller-Hartmann (28). This Hartree term is already 
contained in the LDA. 

In this Section, we will review recent approaches which start from the ah ini- 
tio Hamiltonian [2] and merge conventional bandstructure approaches (Section 
2) with DMFT (Section 3) in order to reliably account for the local correlations 
induced by the local Coulomb interaction. In Section 4.1, we will discuss how 
conventional LDA calculations are merged with DMFT, presenting two points 
of view: (i) a physically-motivated Hamiltonian formulation in Section 4.1.1, 
which can be done self-consistently as pointed out in Section 4.1.2 and whose 
parameters can be calculated by constrained LDA as discussed in Section 4.1.3; 
and (ii) in Section 4.1.4 a formulation in terms of a spectral density functional 
theory which appreciates the Hohenberg-Kohn (2) extremum principle for the 
energy. Section 4.1.5 is devoted to simplified LDA-I-DMFT implementations for 
transition metal oxides. In Section 4.2, the Hartree-|-DMFT and the Hartrcc- 
Fock+DMFT approach are introduced. The former is very natural from the 
DMFT point of view concerning non-local interactions. Section 4.2 also paves 
the way for the more profound alternative to LDA-I-DMFT which is a com- 
bination of the GW approximation and DMFT, finally addressed in Section 
4.3. 



4.1 LDA+DMFT 
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4.1.1 Hamiltonian formulation. In Section 2.3, we deduced the LDA 
Hamiltonian [11] from the ab initio Hamiltonian [2]. To account for the impor- 
tant local correlations, we will now supplement this LDA Hamiltonian [11] by 
the local Coulomb interaction. The resulting multi-band many-body problem 
will then be solved by DMFT. A predecessor of this LDA-I-DMFT approach 
is the so-called LDA-I-U method of Anisimov et al. (19), which constructs the 
same many-body problem but then solves it by (symmetry-broken) Hartree- 
Fock, i.e. without taking into account electronic correlations if these are de- 
fined as being beyond Hartree-Fock. As an approximate solution of the DMFT 
self-consistency equations we will rediscover this LDA-I-U scheme in Section 



In principle, one could use the most general local Coulomb interaction as in 
Hamiltonian [15] of Section 3 and calculate the Coulomb matrix elements by 



However, this approach is not practical. First of all, it is often sufficient to 
take into account only the most important matrix elements of the Coulomb 
interaction, including only two distinct orbital indices and restricting the 
interaction to the localised d or f orbitals only. Another aspect which 
makes Eq. [33] impractical is screening: If we add an extra electron to 
the d or / orbitals at our local site i, the sp electrons will redistribute 
and, therefore, reduce the energy required to add this extra electron. Since 
non-local Coulomb interactions as well as local interactions between d and sp 
electrons are responsible for this screening process, a Hamiltonian with only 
a local Coulomb interaction like Eq. [15] does not describe these processes 
anymore. Therefore, the effective local Coulomb interaction has to be smaller 
(screened) than the overlap integral Eq. [33]. Presently, two methods are 
employed for the calculation of these screened local Coulomb interactions: 
the constrained LDA method in the context of LDA-I-DMFT and the 
random phase approximation (RPA) within the GW-I-DMFT approach. This 
will be discussed in more detail in Section 4.1.3 and in Section 4.3, respectively. 

Set-up of the LDA+DMFT Hamiltonian- Based on these considerations, the 
starting point of the LDA-I-DMFT approach is usually the following Hamilto- 



5.1. 




(33) 
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Figure 10. Pictogram of the intra- and inter-orbital repulsion U and V and the Hund exchange J 

for two orbitals in Hamiltonian [34], 



nian, where i^LDA denotes the LDA Hamiltonian [11]: 



la 



l&Lu lGLu7^m£Lu 

J ' ^ t ' t t 

"■^ ^ila^ila'^ima'^ima ~ ^ila^ila^ima^ima ~ ^j«o-(34) 

ia ia ia 

l£Luj^m£Lu l£Luj^m£Lu l£Lu 



Here, the prime on the sum indicates that every term of the sum is counted 
only once; a =| (j) for a =| (|); n^^ = c^/c^; ; S^a' denotes the Kronecker 
symbol. Moreover, / G Ljj indicates that the interaction is only taken into 
account on the subset Lu of the 'U-interacting' orbitals in the unit cell i. For 
example, Lu might denote the set of d (or the t2g) orbitals of a transition 
metal or the 4/ orbitals of rare earth elements. 

The Coulomb interaction terms consist of the intra-orbital Coulomb repul- 
sion U, the inter-orbital Coulomb repulsion V, the Hund exchange term J 
and a pair-hopping term J, see Fig. 10. The last term of Eq. [34] describes 
the contribution of U, V and J already included in Hlda, which has to be 
subtracted to avoid a double counting. This last term led to some criticism. It 
arises, however, very naturally from the constrained LDA calculation, as we 
will discuss below. If the Coulomb interaction is the same for all orbitals, the 
double counting term corresponds to a simple shift of the chemical potential 
and has no effect. This is the case if, for example, only d orbitals are taken 
into account, see Section 4.1.5. 

In Eq. [34], the Hund exchange has been separated into its z-z and the spin- 
flip components which altogether yield a SU(2) symmetric contribution of the 
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form 



Hj = -2JY, 



(35) 



for the coupling of the spins SjjTi between different orbitals l^vn at site i. For 
real-valued wave functions 

(r) , it immediately follows from Eq. [33] that the 
pair-hopping amplitude equals that of the exchange interaction: J = J. Note, 
that purely real wave functions can be chosen if time reversal symmetry is 
not broken, for example by an external magnetic field. The pair-hopping term 
has not yet been included in LDA+DMFT calculations because one commonly 
assumes that configurations are rare in which one orbital is doubly occupied 
while another is empty if the Hund exchange and the hopping terms are in- 
cluded. For quantum Monte-Carlo (QMC) simulations, the spin-flip term of 
the Hund exchange poses a (sign) problem (149) and has, hence, not been 
included in LDA+DMFT calculations with QMC as an impurity solver so 
far. Different ways to overcome this sign-problem are discussed in Section 5.5. 
These improved methods, such as the continuous time QMC, may hence allow 
to do such LDA+DMFT calculations in the future. 

For degenerate orbitals, e.g. the t2g orbitals in a cubic symmetry, Hamilto- 
nian [34] has to be invariant under orbital rotations, necessitating V = U — 2J. 
Then, with M interacting orbitals, the average Coulomb interaction is 

- ^ U + {M-1){U- 2J) + [M-l)(U- 3J) 

2M-1 ■ ^ ^ 

Even if degeneracy is only fulfilled approximately one often employs V = 
U — 2J as an approximation, which is often quite good and allows us to set up 
Hamiltonian [34] with only three parameters: the average Coulomb repulsion 
U, the exchange interaction J and the double counting correction Ae. All these 
parameters can be calculated by constrained LDA, see Section 4.1.3. 

After setting up the multi-band many-body Hamiltonian [34], we have to 
solve it. This can be done, without genuine electronic correlations, in the 
LDA+U approach of Anisimov et al. (19) or more sophisticatedly by DMFT, 
as in the pioneering work of Anisimov et al. (31) and Lichtenstein and Kat- 
snelson (32). Within DMFT, Hamiltonian [34] is mapped onto an auxiliary 
Anderson impurity problem with the same Coulomb interaction and an ef- 
fective non-interacting Green function {Q^)~^ = (G~^ + 5]), see Section 3. 
This Anderson impurity model has to be solved self-consistency together with 
the k-integrated Dyson equation [32], relating the local Green function G and 
self energy 5]. The LDA part Hlda of Hamiltonian [34] is only entering in 
this Dyson equation [32]. Most conveniently, Hi^da is formulated in k-space 
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with matrix elements as in Eq. [13], and the double counting correction Ae is 
included in this term: 

^LDA = ^LDA- J2 ^^f^il<r= E ^fmH^) 4^^^m (37) 
icrlELu 'k.alm 

with eP°^(k) = eF^^(k) - 5im^e "il e Lu (for all [/-interacting orbitals) 

and e["^ (k) = e]^^^ (k) yi^Lu (for all non-interacting orbitals). Employing 
e^°^(k) of Eq. [37] for the k-intcgratcd Dyson equation [32], the LDA+DMFT 
approach consists of the following steps, see flow diagram Fig. 11: 

(i) We start with a conventional LDA calculation, as in the flow diagram Fig. 
4. 

(ii) With the given density we can calculate the hopping matrix elements tn 
and from these, via a Fourier transformation, the bandstructure ef^^- 

(iii) In a constrained LDA calculation we determine U, J and Ae. 

(iv) This defines the Hamiltonian [34] which we solve by DMFT, see flow di- 
agram Fig. 8. Thereby, ef^^ enters in the Dyson equation; U and J as 
interaction parameters in the Anderson impurity model. 

(v) Prom the DMFT self energy S(a;), we can calculate a new electron density. 
With this, we can go back to step 2. to calculate a new LDA Hamiltonian, 
if we do a self-consistent LDA-I-DMFT calculation as will be discussed now. 



4.1.2 Self-consistent LDA+DMFT calculations. In the first four steps 
of the LDA+DMFT scheme above, the LDA bandstructure calculation and 
the inclusion of electronic correlations by DMFT are performed sequentially. 
In general, the DMFT solution will result in changes of the electron density 
p{r). see Eq. [38] below. The new p{r) then leads to a new LDA-Hamiltonian 
-f^LDA [11] since both the Hartree term and the exchange correlation term of 
Hlda depend on p(r). Also, the Coulomb interaction U changes and needs 
to be determined by a new constrained LDA calculation, resulting in a new 
many-body problem. 

This gives an iteration scheme as in the flow diagram Fig. 11, which without 
Coulomb interaction reduces to the self-consistent solution of the Kohn-Sham 
(60) equations since 5] (a;) = 0. A similar self-consistency scheme was employed 
by Savrasov and Kotliar (150) in their calculation of Pu, albeit without self- 
consistency for U, and without J and Ae. The quantitative difference between 
non-self-consistent and self-consistent LDA-I-DMFT calculations depends on 
the change of the number of electrons in the different bands after the DMFT 
calculation. This change, of course, depends on the problem at hand. For 
example for the Ce calculation presented in Section 6.1.1, the change of the 
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Do a LDA calculation as in flow diagram Fig. 4, Section 2.3: =J> p(r). 




Calculate the tujmS (Eq. [12]) of the LDA Hamiltonian [11]: 

ki jm= f d^ripii{r) + Vion (r) + 
J L 2me 

J <fr'p{v')V^{v-r') + ^^^^]<^,.^(r), 

and the Fourier transform of tuj^, i.e. the LDA band structure 
€LDA(k) (Eq. [14]). 


Do a constrained LDA calculation, see Section 4.1.3: 

Eq.[42) Eq.[50) Eq.[44) 

p(r) — > U; p(r) — > J; p{r) — > Ae. 


Solve the so-defined many-body problem Eq. [34] by DMFT as in 
flow diagram Fig. 8 with the parameters e'"'-'^(k), U, J and Ae. 


Calculate a new density pnew via 

^/ N Eq.[40] Eqs.[39],[38] 
2^{LU) > ^iljm Pnew 


Iterate with p{r) = pnew(r) until convergence, i.e. ||p(r) — Pnew(r)||<e 



Figure 11. Flow diagram of the LDA+DMFT algorithm. 



orbital occupation (from LDA to LDA+DMFT) was very minor in the vicinity 
of a-7 transition but more significant at lower volumes. Hence, a self-consistent 
LDA+DMFT calculation does not appear to be necessary in the vicinity of the 
Ce a-7 transition. Indeed, most LDA+DMFT calculations reported so far are 
not self-consistent. A notable exception is the work of Savrasov and Kotliar 
(150) on Pu already mentioned. 

We still have to describe how to calculate the new electron density p(r) after 
the DMFT calculation. We can relate p(r) to the DMFT self energy E(a;) as 
follows: 

Expressing p(r) via an expectation value of the field operator ^~^{r,a) and 
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doing the basis transformation to the orbital basis set (fiii{r) (Eq. [10]), we 
obtain 

p(r) = 5^(*+(r,a)^(r,a)) = ^ui^Mmirm! &U) ■ (38) 

c cril jm 

Here, the expectation value can be calculated from the equal-time Green func- 
tion 

Gtijmir = 0+) = -5ij5im + {c^i^c'^J. (39) 

This Green function is, in turn, the Fourier transformation with respect to 
space and (imaginary) time of 

G^iiu,) = [iuA + Ml - eL°^(k) - E(ia;,)] (40) 

(bold symbols denote matrices in the orbital index). In principle, the fre- 
quency dependent DMFT self energy 'S{iu}^) and the k-dependent LDA e(k) 
can be spin-dependent. In the present presentation however, such a possible 
spin-dependence is not taken into account. We also have not discussed the 
peculiarities of the LDA basis set; rather we assume that the LDA Hamilto- 
nian has been calculated in an appropriate basis set with sufficiently localised 
orbitals. For the self-consistent LDA-I-DMFT calculation, this basis set should 
also be optimised again with the new LDA-I-DMFT density p{r). 

A remark: Since the solution of the DMFT Anderson impurity model is 
often the computationally most expensive task, it might be a good idea to 
do only one DMFT iteration for every iteration loop in the flow diagram Fig. 
11. When converged, this will produce the same result as a scheme in which, 
for every new p(r), the DMFT equations are iterated until convergence. One 
might also omit the initial LDA loop. However, we would expect that the 
overall numerical effort grows. 

4.1.3 Constrained LDA calculations of the interaction parameters. The 

task of calculating U, J and Ae in Hamiltonian [34] is not at all trivial as 
pointed out on p. 26 and requires additional approximations. Such an ap- 
proximation, which works in practice and gives very reasonable values for 
the screened Coulomb interaction, is the constrained LDA method by Ded- 
erichs et al. (151), McMahan et al. (152) and Gunnarsson et al. (153). Hence, 
constrained LDA allowed us to perform parameter-free, i.e. truly ah initio, 
LDA-I-DMFT calculations for transition metal oxides and / electron systems, 
yielding results in good agreement with experiment (154; 155). 
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Figure 12. Left: In the constrained LDA calculation the interacting d or / electrons on one lattice 
site are kinetically decoupled from the rest of the system, i.e. they cannot hop to other lattice sites. 

By contrast, the non-interacting electrons can still hop to other sites (indicated by dashed lines) 
and screen the d or / electrons. Right: This allows us to change the number of interacting d or f 
electrons on the decoupled site and to calculate the corresponding LDA energies E{n) [circles]; the 

interacting Hamiltonian Eq. [34] reproduces these constrained LDA energies if U and Ae are set 
properly; the dashed line sketches the behaviour of Hi^oA defined in Eq. [11]. 

The basic idea of constrained LDA is to do LDA calculations for a slightly 

modified problem, i.e. a problem where the interacting d ox f electrons of one 
site are kinetically decoupled from the rest of the system. This is achieved 
by setting the hopping matrix elements between these localised orbitals I on 
site i and all other orbitals to zero: tuj^ = \/j,m, see Fig. 12. This allows 
us to change the number of interacting electrons on this decoupled site; they 
cannot hop away. At the same time, screening effects of the other electrons are 
taken into account, since they redistribute if the number of d or / electrons is 
changed on the decoupled site. 

Usually, we know how many d or f electrons Ud to expect, for example, 
from the formal oxidation, and typically Ud is close to an integer number. 
Then, we can do constrained LDA calculations for Ud — 1, Ud, and + 1 
electrons on the decoupled site within the interacting orbitals (/ G Ljj), This 
leads to three corresponding total energies, see Fig. 12. The LDA Hamiltonian 
Hlba which was calculated at a fixed density p{r), see Section 2.3, would 
predict a linear change of the LDA energy with the number of interacting 
electrons, i.e. E{nd) = Eq + e^^^ ra^. This linear behaviour does not take into 
account the Coulomb interaction U which requires a higher energy cost to add 
the (n + l)th electron than to add the nth electron. This effect leads to the 
curvature in Fig. 12 and is taken into account in the Hamiltonian [34] which 
yields for J = 0: 



Eirid) = Eq + \ Unaina - 1) + (e^°^ + Ae)nd. 



(41) 
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Here, (^^^ is the energy contribution of H^ba in Eq. [34] for d or f orbitals on 
the decoupled site. If m denotes the (degenerate) d or f orbital this contribu- 
tion is given by the diagonal hopping matrix element (Eq. [11]): e^^^ = timim 
(if so, orbitally averaged). Note that all hopping elements timji to other or- 
bitals and sites are set zero in the truncated problem of the constrained LDA. 
Since the electrons are kinetically decoupled from the rest of the system, we 
can calculate the interaction energy in Eq. [41] directly from the number of 
electrons.^ Note that part of Ae just arises to cancel the Coulomb contribution 
(1/2) diUnairid - l)]/dnd = U{nd - 1/2). 

Knowing how the energy E{n) changes when going from — 1 to and 
n^ + l electrons allows us to determine U and, Ae since there are two param- 
eters and two energy differences. In the actual calculation, we obtain these 
energy differences from the constrained LDA eigenlevels of the decoupled d or 
/ orbitals, using the relation ddirid) — — dna^ (156j 157). From this we directly 
determine 

U = E{nd + 1) + E{nd - 1) - 2E{nd) 

~ ^d{nd +\)- ^d{nd - ^)- (42) 

Ae is obtained from the energy difference to add a single electron. Eq. [41] 
yields 

E{nd) - E{nd - 1) = U{nd - 1) + (e^^^ + Ae) (43) 

so that 

e^°^ + Ae = E{nd) - E{nd - 1) - U{nd - 1) 

Medina -^)-U{nd-l). (44) 

This shift of the d or f eigenlevel can be taken into account by replacing 
gLDA ^ timim in i^LDA [34] by the term ed{nd - ^) - U{nd - I) from the 
constrained LDA calculation. This replacement thus includes Ae in Elba- 

With this procedure, we guarantee that Hamiltonian [34] correctly repro- 
duces the constrained LDA energies for the truncated problem where hopping 
from and to the localised dor f orbitals onto one site i is forbidden. It is slightly 
different from the approach by Anisimov et al. (19) who directly related the 



^Operators Ylm'"'^ constant all the time and can be replaced by this constant, i.e. n^. 
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double counting correction to U: 

^e = U{na-\). (45) 

A third path has been proposed by Lichtenstein et al. (158; 36), identifying 
an orbital-dependent term 

Ae^ = ^^E^^(a; = 0). (46) 

a 

A cautionary remark: While the total LDA spectrum is rather insensitive to 
the choice of the basis, the constrained LDA calculations of U depend more 
strongly on the shape of the orbitals which are considered to be interacting, 
e.g. the d or f orbitals. For example, in the case of LaTiOs at a Wigner- 
Seitz radius of 2.37 a.u. (atomic units) for Ti, a LMTO-ASA calculation by 
Nekrasov et al. (159) using the TB-LMTO-ASA code (160) yielded U = 4.2 
eV while a ASA-LMTO calculation within the orthogonal representation by 
Solovycv et al. (161) gave U = 3.2 eV. Thus, an appropriate basis is mandatory, 
and, even then, a significant uncertainty in U remains. For conventional LDA 
bandstructure calculations, the choice of basis is not so crucial since with a 
large enough basis set the same bandstructure e^^^(k) will be obtained. In 
this respect, LDA+DMFT is more challenging since it now becomes important 
for which orbitals the Coulomb interaction is taken into account. 

Similarly as U, the Hund exchange J can be obtained from constrained LDA 
calculations. Instead of the number of electrons, however, the spin-polarisation 
has to be changed. Since the electrons within the interacting orbitals on our 
decoupled lattice site cannot hop, we can change the number of n^-f and nai 
electrons, while keeping = n^i + Udi fixed. Then the Coulomb repulsion 
U will not be affected, but nonetheless ii^(nrf|,nrfj = — n^^) will change 
because of the Hund rule exchange. Instead of the change of the total energy 
with n^-f, we can, as in Eq. [42], use the now spin-dependent energy levels 

^dai^dajnaa) = ^^^^^7^^- From this Anisimov et al. (19) determine J as 

/ IT'd ^ ™d Ix / 1 Ix 

-/cLDA = edTlM = Y + 2' M = Y - 2'' ~ ^''i^^'^T = Y^2''^^^^ Y~ T' 

(47) 

Taking such a value of J, Anisimov et al. (19) assume that the Hund exchange 
enters via a term 



Hj = -J ^ niaflma- 



(48) 
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Nature and the LDA+DMFT Hamiltonian [34] are, however, different. As 
depicted in Fig. 10, the Hund exchange also influences the difference between 
inter-orbital and intra-orbital Coulomb interaction, given hy V = U — 2J in 
Hamiltonian [34]. This goes beyond Eq. [48] which only gives a reduction of 
the inter-band interaction for two spin aligned electrons. 

For example, let us consider the case of a V^"*" ion in a cubic crystal field. 
Then we have two d electrons (n^ = 2) in three degenerate t2g orbitals. These 
two electrons can be spin aligned or not. For Hamiltonian [48], the energy 
difference between these two configurations is simply 

S(ndT = 2, ridi = 0) - E{nd^ = nai = 1) = J. (49) 

In contrast, the difference also depends on the orbital configuration for the 
LDA-hDMFT Hamiltonian (Eq. [34]). For the configuration (n^-f = 2, nai = 0), 
the two electrons have to be in different orbitals and their energy is V — J = 
U — 3 J. For the configuration (n^i = l,ndj^ = 1) aiid without the spin-flip 
term, we have the possibility that the two electrons are in the same orbital 
with E{n^^ = rifii = 1) = U or in different orbitals with E(nd] = n^i^ = 1) = 
V = U — 2 J. Averaging over the (unknown) orbital configurations, this yields 
an average value E{nci] = n^i = !) = [/ — 4/3J. Hence, we have 

E(nd^ = 2, nai = 0) - E(nd^ = nai = 1) = (50) 

for Hamiltonian [34] if the spin-fiip and pair-hopping term are not taken into 
account. 

For example, the constrained LDA calculation by Solovyev et al. (161) yields 
a value of JcLDA = 0.93 eV for V'^^. To reproduce this result with Hamiltonian 
Eq. [34], the left hand side of Eq. [50] has to be 0.93 eV and, therefore, the 
parameter J in Eq. [50] has to be J = 3/5 0.93 eV= 0.56 eV. Such a reduced 
value of the Hund exchange coupling is also more reasonable when we compare 
it with the atomic value of J = 0.689 eV which is precisely measurable from 
the atomic spectrum^, see Table 1. In particular, note that it is implausible 



^How to extract the Racah parameters from the atomic spectra is discussed e.g. by Griffith (162). 
He arrives at the following fit formulas for the Racah parameters B and C: 

B=^ [E(3P) - eCf)] + ^ [fi(iG) - E(iD)] (51) 

C = ^ [21E('^F) - SiEi^P) + 'IOECG) + 37ECd) + 28E{'-S)] (52) 

where E(X) denotes the energy of the atomic state X. Loolcing up tlicsc atomic en- 
ergies for the vanadium 3d^ configurations at tlic NIST Atomic Spectra data base 
(http://physics.nist.gov/cgi-bin/AtData/main_asd, which stem in turn from the overview ar- 
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constrained LDA Jclda 



J — 5'/cLDA 



atomic value 



0.9 eV 



0.6 eV 



0.7 eV 



Table 1. Hund exchange for a V^'*' ion as calculated by constrained LDA, corrected value for Eq. [34] and 
the atomic value from Suga and Corliss (163). 



that the Hund exchange is larger than the atomic value since we expect that 
screening effects reduce J in the crystal. 

Similarly, the constrained LDA JcLDA has to be corrected by J = ^ JcLDA 
for n^i = 3 before it is included in Hamiltonian Eq. [34]. These two factors 
determine all correction factors for the exchange coupling J between t2g or- 
bitals because of the unimportance of J for n,^ = 1 and because particle hole 
symmetry gives the same correction factors for = 4 and 5 as for = 2 and 
1, respectively. 

4.1.4 Spectral density functional theory formulation. Our previous 
derivation of the LDA+DMFT method was physically motivated. We started 
from the assumption that the Kohn-Sham equations, i.e. the LDA part, yield 
the correct results for the weakly correlated s or p bands, while the DMFT- 
part takes into account the local Coulomb interactions of the strongly corre- 
lated d or f bands. Let us now consider an alternative way to formulate the 
LDA+DMFT method as a spectral density functional theory, which goes back 
to Chitra and Kotliar (164) and Savrasov et al. (165; 150), for more infor- 
mation also see Ref. (37) and Ref. (36) and particularly the recent review by 
Kotliar et al. (40) . The basic idea is to replace the Hohenberg-Kohn (2) DFT 
functional E[p\ by an energy functional E[p, G] which depends on the electron 
density p(r) and the local Green function Giii„i{uj) of the interacting orbitals 
(m, / G Lif). Since the local Green function is related to the k-integrated spec- 
trum Savrasov et al. (165; 150) coined the name spectral density functional 
theory. 

As a starting point let us take the standard Luttinger-Ward (166) many- 
body functional, which in compact form (suppressing site i, orbital I, spin a 
and Matsubara frequency Ui, indices) reads 



Here Tr denotes the trace (-^ '^uila) ^^'^ ^[^] ^^^^ of two-particle irre- 
ducible diagrams with lines G and without external legs. From the stationary 



n[j:,G] = - TrSG - Trln((G°)-^ - S). 



(53) 



tide by Sugar and Corliss (163)), we obtain B = 0.1085 eV and C = 0.4181 eV so that the Hund 
exchange coupling is J = |-B + C = 0.689 eV for the free atom. 
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conditions 

M = (54) 

we obtain the diagrammatic representation of the self energy: 



6G 

The second stationary condition 



(55) 



= (56) 



(5S 

yields the usual Dyson equation, i.e. Eq. [32] without k integration: 

G = [(G°)-^ - S]-^ (57) 

Such a construction is then a 'conserving approximation' in the sense of Baym 
and Kadanoff (167; 168). On the basis of this functional theory, albeit with 

a purely local Green function and self energy, Janis (169; 170) formulated a 
generalised coherent potential approximation (CPA), and Janis and Vollhardt 
(77) the DMFT self-consistency approach. 

On the other hand, DFT was also formulated as an effective action by 
Fukuda et al. (171; 172). Now, the electron density p{r) plays the role of the 
Green function G in Eq. [53], and the Kohn-Sham potential Vks(r) replaces 
the self energy. In terms of these variables, the DFT functional reads 

i^DFT[FKS,p] = J dVMon(r)p(r) + ^ J d^rdV^^^+E^M 

- J d^r VKs{r)p{r) - TVln[^a;, + h''/2meA - FKs(r)].(58) 

Here, the first line corresponds to the many-body functional (j)[G] which is now 
given by the ionic potential, the Hartree and the exchange correlation term; 
Tr denotes again the trace (now w.r.t. the operator following and any basis, 
e.g. the |r) basis) as well as the sum over Matsubara frequencies as before. 

Minimisation of the functional w.r.t. p, i.e. (5riDFT[p5 ^s]/^P = 0, repro- 
duces the correct form the Kohn-Sham potential in the Kohn-Sham equation 
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[6], 

VKs(r) = Fion(r) + j d\' Fee(r-r')p(r') + (59) 

in analogy to Eq. [55] which gives the correct form of the self energy. The 
second minimisation, ^r^DFxip, Vksl/f^^KS = 0, gives the Kohn-Sham equation 
[6] itself if 14cs(r) is replaced by the right hand side of Eq. [59], which of 
course is guaranteed if self-consistency is obtained. More precisely, it yields 
p(r) directly, expressed via the Kohn-Sham equation: 

p{r) = 1^ (r I [i^. + /iV2meA - Fks] | r) . (60) 

V 

In the case of LDA-I-DMFT, we have to deal with a mixed representation in 
terms of the density p(r) (LDA) and the local Green function of the interacting 
orbitals G(a;) (DMFT). Savrasov et al. Ref. (165; 150) constructed such a 
functional: 

J7[5],G,Vks,p]= j dVMo„(r)/>(r) + i j Srdh' ^^^^ + E^,[p\ 
-$[G]-$dc[G] 

VKs(r)p(r) - TrS(za;^)G(za;^) 
-Trln[ia;^+;iV2"^eA-VKS-S(iu;^)i^^*Kr')^im(r)](61) 



Here, Ti{iu}i,)ijn^*i{r')ipiyn{r) is the transformation of the local orbital self en- 
ergy to the real space formulation in terms of r and r', and $dc[G] is an 
additional double counting correction term. 

In the spirit of Hohenberg and Kohn (2), we are looking for the density 
p{v) and the local Green function G(u;) at which the functional [61] takes its 
minimum. (The minimisation w.r.t. S and Vks gives two additional equations, 
defining SI and Vks-) The minimisation procedure yields the exact ground state 
energy, density and local Green function if £?xc[p] and $[G(a;)] are themselves 
the exact ones. However, we do not know these exact E^c\p\ and $[G(6<;)] 
and have to use approximations. LDA-I-DMFT is one such approximation 
where we replace £'xc[p] by the LDA exchange-correlation energy E^^^[p\ 
and $[G(a;)] by all local two-particle irreducible diagrams, i.e. the DMFT di- 
agrams $DMFT[G(a;)]. The double counting term $dc[G] will have one of the 
forms discussed in Section 4.1.3. 
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With this LDA+DMFT approximation of the spectral density functional 
[61], we recover the LDA+DMFT equations formulated earlier in this Section: 

• The DMFT self energy 5] (a;) minus the double counting correction term is 
obtained from 

• From 

iOP^OWl^O (63) 

and expressing —f^jlnif,^ + Vks in the basis <^ii{y\ i.e. replacing it by the 
orbital matrix (the LDA bandstructure) e^^'^(k) defined via Eq. [14], we 
obtain the Dyson equation: 

G(k, iu,) = [iu, - eLi^^(k) - T.{iu,)]-^ (64) 

• Similarly, 

5J^[5],G,Fks,p] 



= (65) 



yields the Kohn-Sham equations, or more precisely, p(r) expressed through 
the Kohn-Sham equations as in Eq. [60], as in Eq. [60] but now with the 
self energy included: 



p{r) = -^^{^\ [i^u + h^/2meA - Vks - T,{iLOi,)imiPii{r')ipim{r)] 



-1 



(66) 

This equation can also be formulated in the orbital basis fii{r) instead of 
the real space basis |r). 
• Finally, the correct form of the Kohn-Sham potential follows from 

dm,G,VKS,p] ^ Q .g^. 
6p 

as before. 

Let us not forget to mention a conceptional deficiency of spectral density 
functional theory (or LDA+DMFT) compared with standard DFT: What we 
define as local and interacting in DMFT is basis dependent. Different basis sets 
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will hence produce different results, unless more and more orbitals are consid- 
ered as interacting in DMFT and more and more non-local correlations are 
taken into account. This poses the so-far open question: What is the optimal 
basis set for LDA-FDMFT? 

4.1.5 Simplifications for transition metal oxides with well separated eg 
and t2g bands. Many transition metal oxides are cubic perovskites, often with 
a slight distortion of the ideal cubic crystal structure. In these systems, the 
cubic crystal field of the oxygen ions splits the d-orbitals into three degenerate 
t2g and two degenerate Cg orbitals. This splitting is often so strong that the t2g 
or Cg bands at the Fermi energy are rather well separated from all other bands. 
For the low energy physics, it is hence sufficient to take only those bands into 
account which cross the Fermi energy, e.g. the three t2g bands. Note that these 
effective bands at the Fermi energy are not purely of d character. These are 
the effective bands corresponding to the LDA eigenenergies. Hence, they also 
contain contributions of the other orbitals, in particular the oxygen p orbitals 
due to hybridisation effects. 

From a conventional LDA bandstructure calculation, these effective low- 
energy bands can be accurately determined via the NMTO downfolding ap- 
proach of Andersen et al. (46; 47) or the projection to Wannicr orbitals (49). 

If the transition metal oxide is cubic or only slightly distorted, a further 
simplification can be employed for the (almost) degenerate bands at the Fermi 
energy: Instead of the full LDA Hamiltonian Hi,ba in Eq. [34] only the total 
DOS of the low-lying bands needs to be taken into account within DMFT: 
Without symmetry breaking, the Green function and the self energy of these 
bands remain degenerate, i.e. Gim{^) = G{i^)Sim and S;„(ci;) = 'E{Lo)Sim for 
l,m El Ljj (i.e. the interacting orbitals at the Fermi energy). Downfolding to a 
basis with these degenerate Lu bands results in an effective Hamiltonian ^ 
or a corresponding e^^^(k). From this reduced Hamiltonian, the diagonal 
element of the Green function is calculated via 



Due to the diagonal structure of the self energy, the degenerate interacting 
Green function can be expressed via the non- interacting Green function G^{uj): 



Thus, it is possible to use the Hilbert transformation of the unperturbed LDA- 
calculated density of states (DOS) A'''^(e), i.e. Eq. [69], instead of the Dyson 




(68) 




(69) 
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equation (Eq. [32]) with the fuh Hamilton matrix e^^^(k). This simphfies the 
calculations considerably. With Eq. [69] also some conceptual simplifications 
arise: (i) the subtraction of the double counting correction Ae only results in 
an (unimportant) shift of the chemical potential in Eq. [69] and, thus, the 
exact form of Ae is irrelevant; (ii) the theorem of Miiller-Hartmann (66) of 
a fixed spectral function at the Fermi energy holds within a Fermi liquid; 
(iii) as the number of electrons within the different bands is fixed, a self- 
consistent LDA+DMFT calculation (step 5 in Section 4.1.1) is dispensable. 
If the distortion from the cubic symmetry is only weak or if the off-diagonal 
elements of e^^^(k) are small, a slight modification of Eq. [69] can be used: 

G^M = G^iuj - S(u;)) = / de (70) 

J — ^mm\^) — e 

which now depends on the orbitals, albeit neglecting the orbital admixture. 

It should also be noted that the approximation Eq. [69] is justified only if 
the overlap between the t2g orbitals and the other orbitals is rather weak. 



4.2 Hartree+DMFT and Hartree-Fock+DMFT 

Within the concept of DMFT, a combined Hartrcc+DMFT approach is a 
very natural approximation, since the non-local (density-density) Coulomb 
interaction Uij enters only via the Hartree (173) term in the d — > oo limit 
(28): There are equivalent Uij terms so that Uij has to scale like 1/Zji_^j 

in order that the Hartree energy contribution stays finite. This is exactly the 
analogy to the Weiss (174) mean field theory for spin models, in which the 
summed contribution of the neighbouring sites gives an effective Weiss field 
^Weiss _ Y^^_^. Jij{Sj) so that the proper scaling is Jij ~ l/Zjj^jj. 

Consequently, other contributions scale at least like 1 / y/Zjj^jj. For the Fock 
term, i.e. Fig. 13 (right), this scaling is obvious since there is one Green func- 
tion Gij ~ l/y^Zj~~j connecting i and j.^ Note that Si and Smith (175; 176) 
and Chitra and Kotliar (177; 178) employed a different scaling Uj ~ 1/ y^Zjj^jj 
in the so-called extended DMFT, also see Ref. (179). This scaling results in 
additional diagrams, but requires to treat the Hartree term separately since it 
would otherwise diverge. 

The reader might think of the Hartree+DMFT approach as a special case of 
LDA-I-DMFT, in which the LDA exchange and correlation potential Exc[p] = 



^The general way of such scaling arguments was elaborated in Section 3.1. 
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Figure 13. Left: Hartree diagram for the non-local Coulomb interaction. Right: Exchange diagram 
(Fock term). Both are taken into account within Hartree- Fock. 

0. In our orbital basis (pii{r), the Hartree term contributes as 

^flT' = ^Sr^ = / d'r^Ur) I dVp{r')VUr-r')ip,m{r), (71) 

which is the left diagram of Fig. 13 and gives, Fourier-transformed, a k- 
dependent self energy 5]Hartree^ advantage is that Hartree+DMFT is a 
diagrammatic approach, consisting of the Hartree diagram of Fig. 13 and the 
local DMFT contributions of all irreducible diagrams. 

Because we know the Hartree contribution diagrammatically,^ we are in 
full control of the double counting term of the many-body Hamiltonian Eq. 
[34]. Diagram Fig. 13 yields in terms of the on-site equal-time Green func- 
tion G1h^{t = 0^) the following local correction for the self energy of the 
interacting orbitals / G Lfj: 

^.Hartree = ^ [ ^ f;^. = 0+) + ,^(r = 0+) + (F- J)GL .™(t = 0+) 

al '^eit/ 

(72) 

Altogether, the Hartree-|-DMFT self energy consists of the cj-independent, but 
k-dependent, Hartree contribution and the w-dependent, but k-independent, 
DMFT contribution with the double counting term subtracted: 

vHartree-t-DHFT/-, ,\ v^Hartree , r ^DMFTf ,\ r r vHartree /VQ^ 
^ijlm {^) = ^iljm +<^ij^lm [UJ) - dijdlmT.^^ ■ (73) 

With this self energy, one needs to calculate a new electron density p(r) which 
in turn gives a new '^^fj^'^'^'^^^^'^ and has to be iterated until convergence. 



It also follows from the variational principle as was stressed by Slater (180). 



February 2, 2008 



18:28 Advances in Physics ESCwDMFT'16 



Electronic structure calculations using dynamical mean field theory 43 

The advantages of Hartree+DMFT are that it follows the d ^ oo concept 
and that it is a very clear and diagrammatically controlled scheme. The dis- 
advantage, however, is that by neglecting the LDA exchange and correlation 
term all exchange and correlation contributions need to stem from DMFT 
which only includes the local exchange and correlations. With the success of 
LDA and its superiority over a simple Hartree calculation, this seems not to be 
sufficient. This explains the reluctance, in particular, of the physicists of the 
LDA community, to implement such a scheme. Indeed, Hartree+DMFT cal- 
culations have not been done in the context of realistic calculations; for model 
Hamiltonians see e.g. Ref. (101). A more thorough investigation of how well 
such a relatively simple scheme works for realistic calculations is mandatory. 
Also note that with cluster DMFT calculations, which were briefly discussed in 
Section 3.2, more and more non-local exchange and correlation contributions 
can be taken into account so that the treatment of the remaining Coulomb 
interaction terms by the Hartree approximation is a less severe approximation. 
The calculation of the screened Coulomb interaction parameters, however, re- 
mains challenging. 

If one combined the Hartree-Fock (181) approximation and DMFT, at least 
the non-local exchange (given by the right diagram of Fig. 13) would be ad- 
ditionally included. The problem with the exchange term is that the Kohn- 
Sham equations become much more complicated due to the non-local nature 
of the exchange term: In terms of the first quantised wave functions <Pi{Tc), the 
Hartree-Fock equations read 



^ A + V-on(r) + / d'r' ^-^P(r' 
J Aneo r - r 



2mf 



</2j(r) 



d'r- " . > ,(r)=£,99,(r), (74) 



instead of the Kohn-Sham equation [6]. But because it is not local in r, Eq. 
[74] is mathematically much more complicated to solve. Nonetheless, with 
increasing computational resources exact exchange calculations have become 
possible and were performed by Stadele et al. (182). 

Closely related are also the Hartree-Fock calculations by Schnell et al. 
(52; 53). The authors start with a conventional bandstructure calculation 
which is, however, only used to construct maximally localised Wannier or- 
bitals by the Marzari-Vanderbilt (49) projection approach. With these or- 
bitals Schnell et al. (52; 53) calculate the non- interacting hopping matrix 
elements, i.e. the overlap integral of these Wannier wave functions w.r.t. 

A + Mon(r), and the Coulomb interaction by the overlap integral Eq. 
[33]. This defines a multi-band many-body Hamiltonian for which, because of 



February 2, 2008 18:28 Advances in Physics ESCwDMFT'16 



44 K. Held 




Figure 14. The bare Coulomb interaction (wiggled line) is, within RPA, screened by particle-hole 
excitations (Green function bubbles; double lines), resulting in a screened Coulomb interaction 
(wiggled double line) given by an infinite series of bubble diagrams which can be calculated 

self-consistently as indicated. 



the localised orbitals, the hopping and Coulomb interaction between far apart 
orbitals is neglected. Schnell et al. (52; 53) solve this Hamiltonian within a 
reduced set of Wannier orbitals by Hartree-Fock. But their approach also of- 
fers the advantage that more sophisticated many-body approaches like DMFT 
can be employed in the future. The calculation of the screened Coulomb in- 
teraction remains problematic. Schnell et al. (52; 53) propose to employ the 
Thomas-Fermi theory of screening to this end. 



4.3 GW+DMFT 

A serious disadvantage of both, Hartree-|-DMFT and Hartree-Fock-I-DMFT, 
is that the screening of the Coulomb interaction is not yet included. Most 
important for this screening are the random phase approximation (RPA) bub- 
ble diagrams depicted in Fig. 14. Particle-hole excitations, even at relatively 
high energy, dramatically reduce the effective Coulomb interaction. This is 
taken into account by the infinite series of RPA diagrams, replacing the bare 
Coulomb interaction (wiggled line) by the screened Coulomb interaction (wig- 
gled double line). 

For a many-body physicist, a diagrammatically controlled approach for the 
challenging task of realistic material calculations is desirable. The minimal 
set of diagrams to this end are (i) the Hartree and Fock terms of Fig. 13 
which make up a major contribution of the Coulomb interaction, (ii) the RPA 
diagrams of Fig. 14 for calculating the screened Coulomb interaction, and (iii) 
the local DMFT contribution of all irreducible diagrams for d or f materials 
with strong electronic correlations. The former two terms are contained in the 
so-called GW approach, and Biermann et al. (42) recently proposed to include 
the local dynamics of DMFT in a GW-I-DMFT method, also presenting first 
results for ferromagnetic Ni. 

Let us start by briefly recapitulating the GW equations introduced by Hedin 
(41), for a review see Ref. (57). The GW self energy consists of the Hartree 
part, Eq. [71], and the exchange contribution, the right diagram of Fig. 13. 
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The latter is given by 

^G{r,r';iu + oj')W{r,r';iv'). (75) 

This form of the self energy, Green function G times screened interaction W, 
coined the name of the 'GW approximation. The imaginary unit i in front of 
G stems from the standard definition of the real time or real frequency Green 
function and the rules for evaluating the diagram Fig. 13, see Ref. (94). 

To take the RPA screening into account, Eq. [75] employs the screened 
interaction W, i.e. the double line in Fig. 14. This screened W is given by a 
geometric series, defined iteratively in Fig. 14, which yields 

W{r,r';u;) = J £r'%,{r,r")e-\r" ,r' -,0;). (76) 

Here, Vee(r5i"") is the bare Coulomb interaction (Eq. [2] in Section 2] which 
does not depend on lo, and e~^{r,r';Lo) denotes the inverse of the dielectric 
function 

e(r,r';a;) = Srr' - Fee(r, r')P^^(r, r'; a;). (77) 

Together, Eqs. [76] and [77] have the usual pole-structure of the geometric 
series in which the GW polarisation 

P^'^ir, r'; a;) = -2i [ — G(r, r'; uj + J)G{r, r'; lj') (78) 
J 27r 

is the factor for the bubble consisting of two Green functions. In Eq. [78], 
the prefactor 2 stems from the spin summation, assuming here that G and S 
arc spin- independent in a paramagnetic phase. Because of the screening the 
effective interaction W becomes frequency dependent. 

Since the computational burden of the more complicated GW approximation 
is much higher than that of the simpler LDA, plane- wave basis calculations 
become very costly because of the big basis set. Therefore, the LMTO basis 
of Andersen (45) or a Gaussian basis is preferable. 

To merge GW with DMFT, first the local contribution in the GW equations 
has to be subtracted since this part will later be included within DMFT: 

pGW+DMFT(k,^) = pGW(k,^) - -L / d'k^'''^ +P"(a;), (79) 

Vbz J 

BZ 

V ' 

local contribution 
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where we Fourier-transformed ^^^(r, r'; c<j) of Eq. [78] from real to k-space, 
and also employed an orbital basis, indicated by the bold symbols. To switch 
between the representation in terms of spatial coordinates r and the represen- 
tation in terms of orbitals one has to calculate the overlap integral w.r.t. the 
orbital wave functions (pii{r) or to multiply by (r), respectively. In practice, 
Aryasetiawan et al. (44) propose to calculate the two-particle polarisations and 
interactions within another basis, the optimal product basis, see Ref. (57). 
Similarly, we have 



-iGW+DMFT 



VBZ J 



BZ 

^^.Hartree^]^^ O;) - ^""^^ + 5]°^^^T(k, u). (80) 

Here, also the Hartree self energy and its local part, Eqs. [71] and [72], re- 
spectively, have been included. The local DMFT self energy contribution is 
that of the auxiliary Anderson model or, equivalently, that of the local con- 
tribution of all irreducible diagrams. Note, that the combined contribution 
j^DMFT^l^^ w) - x;Hartree actually the exchange-correlation part of the DMFT 
self energy. 

To calculate the DMFT polarisation, let us start with the local impurity in- 
teraction of the DMFT Hamiltonian Eq. [15], considering only density-density 
type of interaction between two orbitals, i.e. Uimno = ^mn^ioUim in Eq- [15]. 
Within DMFT, this local interaction is screened, yielding a local screened 
interaction 

Wim{i~^) = Uim{uj) - ^ Uii'{u)xil'm'm{i^)Um'm{i^), (81) 
m'l' 

which can be calculated from the local susceptibility 

Xlmnoir) = {T4{r)c^{T)ci{0)nM)- (§2) 

Here, T is the time-ordering operator. 

From the difference between the inverse matrices in the orbital index W 
and U, the DMFT local polarisation can be calculated by an inversion w.r.t. 
the orbital indices: 

pDMFT(^) = V-\lo) - W-\u;). (83) 



This relation follows directly from the definition of the local polarisation anal- 
ogous to Eq. [76], after substituting e{io) of Eq. [77] and solving for P(c<;), with 
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U(ti;) playing the role of Vcc- 

Now, however, the interaction U(ct;) is frequency dependent or (in imaginary 
time) r-dependent. This has to be taken into account in the effective action 
of the Anderson impurity model, where we have now to employ the following 
^[V',V'^(a°)-^U] in Eq. [27]: 

va Im 

13 

Im aa' q 

Altogether, we are now in the position of formulating the GW+DMFT 

scheme, introduced by Biermann et al. (42; 43; 44), sec flow diagram Fig. 15: 
Starting point is a conventional LDA calculation, yielding an electron density 
p(r) and also a LDA Green function 

G^^{u) = [ul + /xl - e^^^{k)]-\ (85) 

This Green function is inserted into Eq. [78] to calculate the GW polarisation 
P'^^(k, w) and, via Eq. [77], the screened interaction W(k, w), Eq. [76]. 

Now three self energies have to be calculated: (i) the GW self energy 
I]^'^(k,a;), i.e. the Fourier transform of Eq. [75]; (ii) the Hartree self energy 
and double counting correction, Eqs. [71] and [72], respectively; and (iii) the 
DMFT self energy. For the latter, we first have to define the non-interacting 
Green function of the auxiliary DMFT impurity problem, i.e. Q^. Moreover, 
we have to unscreen the local Coulomb interaction since the DMFT diagrams 
include the local screening bubble, the contribution of Fig. 14 with i' = f. 
This contribution would be doubly counted if we started DMFT with the fully 
screened W. 

The result of the calculation of the impurity problem will be the impurity 
Green function, which is identical to the local GW+DMFT Green function, 
and the local susceptibility. From these two quantities we can determine the 
local DMFT self energy and polarisation. Altogether, this then allows us to 
calculate the total GW+DMFT self energy and GW+DMFT polarisation. The 
former in turn allows us to calculate the GW+DMFT Green function by the 
Dyson equation [24] in which the non-interacting Green function is simply 
the solution of the ab initio Hamiltonian (Eq. [2]) without electron-electron 
interaction. With this new Green function G, the DMFT self energy 5]DMFT^ 
and the GW+DMFT polarisation pGW+DMFT^ restart our iteration 

loop until convergence. 
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Do LDA calculation, Fig. 4, yielding Gk(w) = [wl+^l-eL°^(k)]-i (Eq. [85]). 




Calculate GW polarisation P°^'(w) = -2i/^G(w + (Eq. [78]). 


If DMFT polarisation pDMFT jg known (after the 1 . iteration) , include it 
pGW+DMFT(j^^^)^pGW^k,a;)--^/ d3fcpGW^k,a;)+pDMPT(w) (Eq. [79]). 


With this polarisation, calculate the screened interaction (Eqs. [76], [77]): 

W(k:^0 = V,,,(k)[l - V,,,(k)P(k:^0]"'- 


Calculate T,^^-"-""' (Eq. [71]) and double counting correction S^;-^*-^"" (Eq. [72]). 


Calculate S°^(r, r'; w) = i / |^G(r, r'; w + J)W{r, r'-w') (Eq. [75)]. 


Calculate the DMFT self energy sDMFT ^^^^^ polarisation pDMFT fono^g. 




From the local Green function G and old self energy S^'^^"'" calculate 
(gO)-i(w) = G-i(a;) + S°^''T(a;) (Eq. [31]); S°^''t=0 in 1. iteration. 


Extract the local screening contributions from W (Eq. [83]): 
U(a;) = - pDMFT(^)]-i_ 


With U and Q'^, solve impurity^problem with effective action (Eq. [84]) 
resulting in G and susceptibiUty x (Eq. [82]). 


From G and x, calculate S^mft^^) ^ {g^Y^{uj) - G-'^{uj) (Eq. [31]), 
pDMFT(^) = U-i(a;) - [U- UxU]-i(w) (Eqs. [81], [83]). 


Combine this to the total GW self energy (Eq. [80]): 


From this and G°, calculate GJJ«^(a;)-i = G^(a;)-i - Sk(a;). 


Iterate with Gk = G^*^^ until convergence, i.e. ||Gk - G^*^^! <e. 



Figure 15. Flow diagram of the GW+DMFT algorithm. 



Certainly, such a fully self-consistent GW+DMFT scheme is a formidable 
task, and Biermann et al. (42) employed a simplified implementation for their 
calculation of Ni: For the DMFT impurity problem, only the local Coulomb 
interaction between d orbitals was included and its frequency dependence was 
neglected W(a;) = W(0). Moreover, only one iteration step has been done, 
calculating the inter-site part of the self energy by GW with the LDA Green 
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function as an input and the intra-site part of the self energy by DMFT (with 
the usual DMFT self-consistency loop, see flow diagram Fig. 8). The GW 
polarisation pGW calculated from the LDA instead of the GW Green 
function. This is, actually, common practice even for conventional GW calcu- 
lations. 

Notwithstanding, Aryasetiawan et al. (183) emphasised the importance of 
the frequency dependence of the Coulomb interaction, even for the effec- 
tive low-energy physics. For a model Hamiltonian, the Hubbard model with 
nearest-neighbour interaction. Sun and Kotliar (184; 185) already performed 
such self-consistent GW-I-DMFT calculation and included further local corre- 
lations through the extended (E-)DMFT (176; 175; 186; 178). Also note, that 
one can embed GW+DMFT (similar to LDA+DMFT) as a Luttinger-Ward 
(166) free energy functional which is, here, a functional of the Green function 
G and the screened Coulomb interaction W; for more details see Ref. (43). 



5 DMFT solvers suitable for material calculations 

The equivalence of the DMFT single-site problem and the Anderson impurity 
problem, noted by Georges and Kotliar (29), also see Jarrell (30), allows us to 
employ a variety of well established techniques to solve the DMFT equations. 
Since the DMFT results depend on this solver, we should be careful to employ 
the respective solvers only in the parameter regime where they are applicable, 
and e.g. not a perturbative approximation in the non-perturbative regime. 
Otherwise, we run the risk of reporting a good agreement with experiment 
which stems from a cancellation of errors. LDA+DMFT is an approximation 
and misses non-local correlations which are certainly important for some ma- 
terials. If we introduce a second error due to the DMFT impurity solver, good 
agreement with experiment might be by chance. Possibly, we then overlook 
important non-local correlations. 

For these reasons, Nekrasov et al. (159) coined the notation DMFT(X) and 
LDA-I-DMFT(X) where X denotes the solver for the auxiliary DMFT impurity 
problem. 

In the following we will review the most common techniques X which are 
promising for material calculations with DMFT. At the beginning, we will 
discuss approximate techniques, starting with the Hartree approximation in 
Section 5.1 which makes the LDA+DMFT approach equivalent to LDA+U. 
This allows us to describe an insulating phase with Hubbard bands, albeit 
only if there is spin or orbital order. In the case of alternating antiferromag- 
netic or orbital order, however, LDA+U fails to yield the correct shape of the 
Hubbard bands since it misses the formation of polaron side bands (25). Alter- 
natives, in particular for the paramagnetic phase without symmetry breaking. 
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are the approximations by Hubbard (15; 187) and extensions which are dis- 
cussed in Section 5.2. In the opposite hmit of weak couphng (weak Coulomb 
interaction), one can do iterated perturbation theory (IPT) which is second 
order perturbation theory in the local Coulomb interaction of the Anderson 
impurity model. IPT will be discussed in Section 5.3. 

A good description of the behaviour at strong coupling is possible by 1/A^ 
type of approaches such as the non-crossing approximation (NCA), which is a 
resolvent perturbation in the hybridisation of the Anderson impurity model. 
These approaches become exact for a large number of degenerate orbitals 
and are good for the insulator but also for the metal if the temperature is 
higher than the Kondo temperature of the impurity model. In Section 5.4, we 
will introduce the NCA. 

In Section 5.5, the concept of numerically exact quantum Monte Carlo 
(QMC) simulations is introduced, which directly solve the Anderson impu- 
rity model, albeit on the imaginary axis. Therefore, the maximum entropy 
method (188) is needed for an analytical continuation to real frequencies if, 
for example, the spectral function is calculated. The numerical effort of QMC 
grows cubically with decreasing temperature so that it is possible to do room 
temperature, or say 100 K, calculations nowadays. Since the effort only grows 
cubically (not exponentially) with temperature calculations at lower temper- 
atures will be possible in the future, although not at really low temperatures. 
To this end, Feldbachcr et al. (189) recently developed a projective QMC 
algorithm for 'T=0' calculations which we will also discuss. 

We will not present details of the exact diagonalisation (ED) and the nu- 
merical renormalisation group (NRG) method of Wilson (190), which have 
been very successfully employed for the one-band Hubbard model: ED by 
Caffarel and Krauth (89), also see Ref. (56; 191), and NRG by Bulla et al. 
(192; 193; 194). Hence, let us briefly state the idea and the limitations for 
multi-orbital calculations here. Both of these methods diagonalise the An- 
derson impurity model. In the case of ED the Anderson impurity model is 
diagonaliced directly for a limited number of non-interaction 'bath' sites. As- 
tonishingly, Potthoff (195) demonstrated that even a single "bath' site gives a 
very good estimate of the critical U value of the Mott transition. In the case 
of NRG, one diagonalises within a restricted energy window which is renor- 
malised to lower and lower energies. Since the number of quantum mechanical 
states grows exponentially with the number of orbitals in the Anderson im- 
purity model, these methods are severely restricted concerning the number of 
orbitals: one orbital is manageable; a sound calculation with two orbitals is 
already almost exceeding the computational limits, but can with great numer- 
ical effort still be done; a reliable treatment of three orbitals will be impossible 
even if computer resources grow considerably. Hence the usefulness of ED and 
NRG for realistic material calculations is very limited and we will not dis- 



February 2, 2008 



18:28 Advances in Physics ESCwDMFT'16 



Electronic structure calculations using dynamical mean field theory 51 

CUSS these methods here, referring the reader to the references above for more 
information. 

Most recently, also the related (dynamical) dynamical matrix renormalisa- 
tion group (DMRG) approach (196; 197; 198) has been employed by different 
groups (199; 200; 201; 202; 203; 204). DMRG is a powerful alternative to ED 
and NRG and might become a standard impurity solver for DMFT in the 
future. For realistic multi-orbital calculations, it however also scales exponen- 
tially with the number of orbitals involved. 

Another method which one might subsume under the d — > oo limit is 
the Gutzwiller approximation. As was shown by Metzner and Vollhardt (27; 
205; 206) and Gebhard (207), this approximate treatment of the Gutzwiller 
(208; 209; 210) wave function becomes exact for c? — > oo. Gebhard (207; 211) 
also showed that the d — >■ oo Gutzwiller wave function is equivalent to 
the Kotliar-Ruckenstein (212) slave Boson approximation at zero tempera- 
ture (also note recent slave-rotor (213) and slave-spin (214) variants). The 
Gutzwiller wave function starts with the Fermi sea onto which the so-called 
Gutzwiller correlator g"^ is applied, where d is the operator of local double 
occupations, i.e. d = J2i ^il^il for the one-band Hubbard model. With 5 as a 
variational parameter, this wave function emulates a central effect of the local 
Coulomb interaction, the reduction of the number of doubly occupied sites. An 
important step for our understanding of the Mott-Hubbard transition was the 
work by Brinkman and Rice (20) who recognised that the Gutzwiller approxi- 
mation describes a Mott-Hubbard metal-insulator transition, indicated by the 
disappearance of the quasiparticle peak. This approximation, stemming from 
the metallic side, however misses a correct description of the Hubbard bands. 
Gebhard (211) generalised the Gutzwiller approximation to the multi-orbital 
situation, allowing subsequently Strack and Vollhardt (215; 216) and Gulasci 
et al. (217) to apply it to model Hamiltonians. More recently, it was employed 
for realistic calculations of ferromagnetic transition metals by Biinemann et 
al. (218; 219; 220), for a summary of these results see Ref. (221). Closely re- 
lated to this are also material calculations with the local ansatz, which has 
been employed for material calculations long before DMFT, see, e.g. Refs. 
(222; 223; 224). 



5.1 Polarised Hartree-Fock (HF) approximation (LDA + U) 

The simplest way to deal with the auxiliary Anderson impurity problem of 
DMFT, Eq. [26], is to treat the interaction in the Hartree-Fock (HF) approx- 
imation. The Hartree diagram has already been shown in Fig. 13 (left) and 
corresponds to the decoupling 
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Hence, correlations given by 

(87) 

are neglected within the Hartree approximation, i.e. within the approximate 
expectation value of the right hand side of Eq. [86]. In principle, we can treat 
the Fock term on a similar footing with a decoupling into non-diagonal ex- 
pectation values of the form (c^^cj^). However, for the local Fock term to 
contribute we need either a (non-diagonal) magnetisation in the xy plane or a 
corresponding kind of orbital ordering. At least for degenerate orbitals, these 
Fock expectation values, or the off-diagonal Green function required in the 
diagram Fig. 13 (right), are zero. Therefore, the local Fock term is usually not 
considered, only the Hartree term. 

If the Coulomb interaction has the form of Eq. [34] with an intra-orbital 
Coulomb interaction U , an inter-orbital interaction V and Hund's exchange 
coupling J, Eq. [86] or equivalently the Green function of Fig. 13 (left) result 
in the Hartree self energy 

^f-^^^^ = U{h!)+ V{nl) + {V-J){nl,). (88) 

This is (up to a constant) the same as Eq. [72] if orbitally averaged and 
expressed in terms of Green functions instead of expectation values (Eq. [39]). 

Doing such a LDA+DMFT(HF) calculation is completely equivalent to the 
LDA-I-U approach (19), also see Ref. (225). For a simple Hartree calculation, 
the DMFT formulation is not necessary and instead of decoupling the in- 
teraction of the Anderson impurity model one can directly apply the Hartree 
decoupling to the LDA-constructed many-body multi-orbital problem Eq. [34] . 
This is how Anisimov et al. (19) arrived at the LDA-I-U scheme. 

Because of the purely local interaction, the Hartree self energy, Eq. [88], is k- 
independent. Moreover, it is static, i.e. u;-independent. What kind of physics 
can we then expect to be described? For the paramagnetic phase with an 
equal orbital occupation, sp^''*'''^'^ is independent of /. Hence, the self energy is 
reduced to a constant shift like a chemical potential, which moreover should 
be canceled by the orbitally-averaged double-counting correction given by Eq. 
[72]. There is no effect at all. 

This changes if the spin or orbital degrees of freedom are ordered (polarised) , 
either in a homogenous way like the Stoner (226) ferromagnet or in a more 
complicated pattern such as the checkerboard alternating Neel (227) state, 
which gives rise to the Slater (96) bands. In these polarised phases, the bands 
of different spin or orbital species split into subbands since the self energy Eq. 
[88] differs. For a large average Coulomb interaction U = [U + (M — 1)V + 
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(M— 1){V — J)]/{2M — 1) (M: number of orbitals) this splitting results in two 
completely separated sets of bands with a typical distance at large U. For an 
integer filling, the lower set of bands will be completely filled and the upper set 
completely empty, resembling the Hubbard bands of a Mott insulator, shown 
schematically in Fig. 2 of Section 1. 

Since the number of double occupations is minimal for this large U phase, 
the symmetry-broken Hartree-Fock solution also correctly produces the total 
energy in the U ^ co limit. The corrections to this limit are of order W'^/U'^ 
and depend on the kind of ordering, leading to the LDA+U prediction for 
the symmetry breaking in a specific material. Together with the insulating 
spectrum, this explains why the LDA+U method has been successfully applied 
for the ab initio calculation of insulators, including their low temperature 
ordering. 

This good description of the insulator might be astonishing if one considers 
that the starting point, the Hartree self energy Eq. [88], is first-order pertur- 
bation theory in the Coulomb interaction, which suggests that Hartree is good 
at weak coupling (Coulomb interaction) only. However, except for almost un- 
correlated systems with U jW <C 1 , the Hartree approximation is rather bad 
for the metallic phase. The reason is that it strongly overestimates the energy 
of the paramagnetic phase. This paramagnetic energy grows ~ IJ since dou- 
ble occupations cannot be avoided without symmetry breaking in the Hartree 
approximation. Therefore, the tendency towards orbital or magnetic ordering 
is heavily overestimated, and LDA-I-U almost automatically yields split-bands 
even if this is not correct at small and intermediate strengths of the Coulomb 
interaction where the system is paramagnetic or has a much smaller magneti- 
sation than the Hartree approximation predicts. 

Nonetheless, a paramagnetic insulator with Hubbard bands with separation 
tj can be described by LDA-I-U if one considers the symmetry breaking as an 
artificial means to produce the correct spectrum and energy. Of course, then 
the magnetisation is incorrect, but also the free energy since the entropy of 
paramagnetic uncoupled spins is much higher than that of ordered spins. These 
obstacles can be overcome by the Hubbard-I approximation which also has the 
advantage that the weights of the Hubbard bands are better reproduced for 
non-integer fillings, see the following Section. 

A cautious remark is also in place if one applies LDA-I-U for studying antifer- 
romagnetic or orbitally ordered phases: While LDA-I-U then yields the correct 
static properties, dynamic properties are completely wrong. This was shown 
by Sangiovanni et al. (25), who particularly pointed out that the DMFT Hub- 
bard bands are very different from the coherent, narrowed Hubbard bands of 
LDA-fU: The DMFT Hubbard bands are incoherent, have spin-polaron side 
peaks and their widths correctly converges to that of the non- interacting bands 
for [/^ oo (25). 
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5.2 Hubbard-I, Hubbard-III and alloy- analogy approximation 

Starting point of the Hubbard (15; 187) approximations is the atomic Hmit 
(yV = 0), in which we can calculate the (purely local) Green function exactly, 
including the full multiplet structure with the spin-flip and pair-hopping term 
of the LDA+DMFT Hamiltonian (Eq. [34]). Due to the itinerancy of the sys- 
tem these atomic levels are, however, broadened. This broadening is included 
in different ways in the Hubbard-I (15) and Hubbard-III approximation (187), 
as well as in the alloy-analogy approximation. 

The Hubbard approximations can be derived by decoupling the higher order 
Green functions in the equations of motion. This procedure is not controlled 
and allows for many different decoupling schemes, i.e. approximations. For 
infinite dimensions, the simplified momentum summations allow for the de- 
coupling of higher order Green functions (228). 

For realistic material calculations, Lichtenstein and Katsnelson (32) for- 
mulated the LDA-I-Hubbard-I, or in our notation LDA-I-DMFT (Hubbard-I), 
approach as one of their LDA++ approaches. For the 7-phase of cerium such 
LDA-hDMFT (Hubbard-I) calculations by McMahan et al. (229) were quite 
successful. These cerium results will be presented in Section 6.1.1. 

Here we will motivate the Hubbard-I approximation physically without 
attentiveness to the original equation of motion derivation by Hubbard 
(15; 187; 14) In the atomic limit, the hybridisation of our DMFT Anderson 
impurity model (Eq. [26]) vanishes and the impurity problem is simply the 
isolated atom for which the Green function is known exactly. Expressed via 
the spectral representation, see Ref. (94) p. 163 ff., it reads 



Here, |^) are the exact many-body eigenvectors with eigcnencrgies 

{En) for A'^ [N + 1) electrons. For the atomic limit we know these eigenenergies 
and vectors so that, without Hund's exchange {J = t) = U = V) and for 
degenerate levels with site energy £, the atomic Green function reads 



{ll\Ci \u){v\c]n\lj) 
U + ^ - En + E^, 



g-/3[i?.-M(^+i)] {I _ g-/3(£;^-£;.+M))_ 1^89) 




(90) 



Thereby, the weights of the poles are given by 



(A^ + l)t;jv+i + (2M - Ar)i;jv 
2ME^lo^iV' 



(91) 
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These are determined in turn by the weight for having N electrons on the 
atom with altogether 2M states 

N\{2M-N)\ ' ^ ^ 

consisting of the Boltzmann weight and a combinatorial factor for the number 
of orbital configurations with electrons in M degenerate orbitals. Prom G'^^, 
we can also calculate the average number of local electrons n^*, for example, 
via the following sum over Matsubara frequencies Ui,: 



2M T ^ G'^* (iw^)e^'^-°^ + 2M. (93) 



Since the non-interacting Green function of the Anderson impurity model is 
known to be ^"(w) = l/{uj + jj, — s) the exact self energy in the atomic limit 
reads: 

S^*(u;) =u + i^-s- [G''\uj)]-\ (94) 

Hubbard (15) now approximated the self energy of the itinerant problem by 
this atomic self energy: 

S(a;) = (95) 

This self energy allows us to calculate the k-dependent Green function via the 
Dyson equation [24] or the local Green function via Eq. [32]: 

G(a;) = [ <fk [ul + /xl - e^^^{\i) - S(a;)l]-\ (96) 
Vbz J 

BZ 

Here, e^^^ denotes the LDA bandstructure minus double counting correction 
as defined in Eq. [37]. For consistency, wc require that the number of electrons 
in the auxiliary atomic problem ra*^* equals the average number of electrons 
in the interacting d or f orbitals, which is calculated from the corresponding 
orbitals of the Green function G(a;) of Eq. [96]. For this requirement to hold, 
we have to adjust //at correspondingly (230; 231; 232; 233). Noting that the 
Wj's sum to one, we also see that S^^C^) ""-^ (2M-1)/(2M) Una which 
is the paramagnetic Hartree-Fock value and correct in the high-temperature 
limit. 

Let us now elucidate the features of the Hubbard-I approximation. In con- 
trast to the atomic limit, the Hubbard-I Green function has 'dynamic' Hub- 
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bard bands with a finite width. Usually T ^ [/, so that only two addends 
contribute to the atomic Green function of Eq. [90]. Then, there will be two 
Hubbard bands which are centred around — /x and U — jj. at large U. If the 
average number of electrons in the interacting orbitals Ud is integer, the lower 
Hubbard band at will be completely filled and the upper Hubbard band 
&t U — ji completely empty. Similar as in the Hartree-Fock approximation of 
Section 5.1, but now in the paramagnetic phase. In contrast to the rigid bands 
of the Hartree approximation, these spectral weights shift, however, when Ud 
is non-integer. This is very physical, since for non-integer rid one can add an 
extra electron without paying more Coulomb interaction energy than for the 
last electron added. 

The major deficit of the Hubbard approximations is that the metal which 
is generally described at non-integer is not a Fermi liquid. This is a conse- 
quence of the construction with the atomic self energy. Hence, the Hubbard-I 
approximation is not an adequate method for non-integer Ud, except for such 
high temperatures at which the Fermi liquid behaviour has been lost. Or if, as 
in the case of manganites, strong scattering destroys the Fermi liquid behav- 
ior, justifying the application of the Hubbard-I approximation for this system 
(234; 235; 236; 237) to some extend. 

There have been several attempts to improve the Hubbard-I approximation, 
starting with Hubbard's own work. The alloy-analogy approximation and the 
Hubbard-III approximation (187) introduce an additional energy-dependent 
scattering correction A(a;) to the self energy which then becomes 

S(a;) = S^*(a;- A(a;)). (97) 

In the alloy-analogy approximation, A(a;) is fixed by requiring that the local 
Green function calculated from S(a;) via Eq. [96] still fulfils the atomic-limit 
Eq. [94] 

G^\uj) = [uj + iiat-e- A(6c;) - (98) 

but with the scattering correction /S.{lj). After substituting S^* from Eq. [94] 
into Eq. [98], we obtain an implicit equation for A(c<;). 

The physical picture behind the alloy-analogy is that the electrons of one 
orbital move in a background of static electrons of the other orbitals. This 
reduces the Hubbard model to the simpler Falicov-Kimball model. With 
this approximation, the Coulomb interactions are reduced to a static scat- 
tering potential which is evaluated within the coherent potential approxi- 
mation (CPA). Note that for disordered (non-interacting) systems the CPA 
(238; 239; 240; 241; 242; 243) becomes exact for the scattering problem in 
d ^ oo (244); for a CPA-inspired derivation of DMFT, see (77; 245; 246). 
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Another way of thinking is in terms of moments of the interacting spec- 
tral function. Then, the Hubbard-I approximation yields the first two mo- 
ments correctly and further-going approximations by Nolting and coworkers 
(247; 248; 249; 250; 251; 252) the first three moments. 

The advantage of this alloy-analogy approximation and also the further go- 
ing Hubbard-III approximation axe that the splitting of the bands now occurs 
at a finite U. Moreover, the widths of the Hubbard-III bands at large U equals 
the width of the non-interacting system in agreement with general arguments. 
In contrast, this width is reduced in the Hubbard-I approximation and even 
goes to zero for a Slater (96) type of alternating ordering in the Hartree ap- 
proximation. In contrast, the DMFT Hubbard bands have the correct width 
W also for staggered ordering (25). 

The equation of motion decoupling scheme by Lacroix (253; 254) for the 
Anderson impurity model is also noteworthy here. It results in a three peak 
structure with a central resonance and two split-off side bands. Within the 
DMFT self-consistency scheme, it was employed by Jeschke and Kotliar (255). 

Let us conclude here that Mott insulating materials can be described, at 
strong coupling, by the Hubbard-I and -HI approaches, the alloy-analogy or 
CPA, and the LDA-I-U method. Important differences concern the magnetic 
ordering tendencies and the shape of the Hubbard bands. The development 
of an ab initio approach which can be more generally used than LDA-I-U 
or LDA-I-DMFT (Hubbard-I) for insulating transition metal oxides and other 
Mott insulators seems to be manageable. 

5.3 Iterated perturbation theory (IPT) and extensions 

The approaches discussed so far (Section 5.1 and 5.2) were only eligible deep 
within the Mott insulating phase. To describe the metallic phase, a natural 
starting point is weak coupling perturbation theory in the Coulomb interac- 
tion U . To first order in [/, this is the Hartree-Fock approximation which is, 
however, better suited for a Mott insulator with magnetic or orbital ordering, 
as was shown in Section 5.1. To second order in we have to deal with the 
diagram shown in Fig. 16. As Green function we can either inset the bare non- 
interaction Green function of the Anderson impurity problem or we can 
use the interacting Green function G of the Anderson impurity model, gen- 
erating an infinite series of diagrams which is self-consistent in the sense 
of Baym and Kadanoff (167; 168). Miiller-Hartmann (66) and Schweitzer and 
CzychoU Ref. (82) followed the latter path, which was extended by Menge and 
Miiller-Hartmann (256) to bubble and ladder summations. However, while 
this self-consistent second order perturbation theory gives to order (U/W)'^ 
the same result as the perturbation theory in the bare Q^, it lacks important 
physical aspects, in particular the formation of Hubbard bands. 
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Figure 16. Second order perturbation theory for the Anderson impurity model with the bare 
Green function (left) and the full Green function G (right). 

These and even the correct quahtative features of the Mott-Hubbard tran- 
sition in the one-band Hubbard model are described by the iterated (second 
order) perturbation theory (IPT) of Georges and Kothar (29), in which the 
bare of the auxihary DMFT impurity model is used. Yosida and Yamada 
(257; 258) and Zlatic (259) realised already that this bare perturbation theory 
for the Anderson impurity model is astonishingly good, showing a three peak 
structure with a qualitatively correct central quasiparticle resonance. This 
motivated Georges and Kotliar (29) to employ IPT as a DMFT solver. Later, 
Zhang et al. (260) recognised that IPT, albeit being a weak-coupling perturba- 
tion theory, also becomes correct in the strong coupling limit {U — > oo). Hence, 
it might not be so astonishing that IPT describes the qualitative features of 
the Mott-Hubbard transition correctly, see Ref. (56). 

Let us now write down the IPT equation. From the Feynman diagram Fig. 
16, we obtain in terms of imaginary time r and Matsubara frequencies u^n'- 

P 

SiPT^za;„) = ^-^^n,U -U' j dre^^-^ G\r)G\T)G\-r) . (99) 



^ V ' 

s<^>{j<^„) 

Here, the first part is the Hartree term, which yields the high-u; behaviour. 
The second part stems from the second order diagram of Fig. 16. 

Away from the, maybe coincidentally, good description for the half-filled one- 
band Hubbard model, IPT is less good at intermediate coupling and fails to 
reproduce the atomic limit. Because of these shortcomings but also because of 
the success of IPT for the half-filled one-band Hubbard model, several groups 
developed extrapolation schemes for the self energy which reproduce the cor- 
rect self energy to second order in U /W and the atomic limit. Edwards and 
Hertz (261; 262) extended the Hubbard approximations from Section 5.2 to 
yield the correct result to order [U/W)'^. This approach gives a Fermi liquid 
phase at small C/, however, there is an unphysical paramagnetic, non-Fermi- 
liquid metallic phase at intermediate U values before the system becomes 
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insulating. 

Kajueter and Kotliar (263; 264) made an ansatz for the self energy of the 
form 



,PT', , 2M-1 ^ AJ:('){uj) 



where S*^^^ is the second order contribution from Eq. [99], analytically contin- 
uated to the real axis. The parameters B and A are determined to give the 
correct self energy in the strong coupling limit (U oo) and the correct 1/uj 
behaviour, respectively. The latter also guarantees that the first two spectral 
moments are correct. However, the quasiparticle peak comes out too small in 
the IPT approximation by Kajueter and Kotliar. Potthoff et al. (265; 266) 
employed the same ansatz but adjusted the parameters to reproduce the cor- 
rect third moment of the spectral function as well. In essence, this scheme is 
an interpolation scheme for the self energy between IPT and the strong cou- 
pling spectral density approach by Nolting and Borigiel (248). The ansatz by 
Potthoff et al. (265; 266) also allows for the treatment of magnetic phases. 

Another self energy interpolation scheme was recently proposed by 
Oudovenko, Savrasov et al. (267; 268), using the Kotliar-Ruckenstein (212) 
slave Boson approach at weak and the Hubbard-I approximation at strong 
coupling. Also the local moment approach of Logan et al. (269; 270; 271) 
fulfils this feature of reproducing the correct weak and strong coupling limit. 

Anisimov et al. (31) used IPT with the Kajueter-Kothar IPT (264) inter- 
polation as the impurity solver in the first LDA+DMFT calculations. As was 
demonstrated by Nekrasov et al. (159) in the context of realistic calculations 
for LaTiOs, IPT yields too small quasiparticle weights and violates the Fermi 
energy pinning, which was proven by Miiller-Hartmann (66) for the one-band 
Hubbard model repeating the Luttinger-Ward (166) arguments and which 
generally holds for multi-band systems if the self energy for all orbitals is 
degenerate. This is, of course, generally not the case in realistic material cal- 
culations. It can be fulfilled if the DMFT calculation works with a restricted 
set of orbitals, e.g. degenerate t2g orbitals for transition metal oxides. 

An alternative weak-coupling perturbative approach is the fluctuation ex- 
change approximation (FLEX) of Bickers and Scalapino (272) which has been 
employed as a DMFT solver in realistic calculations by Lichtenstein and Kat- 
snelson (32), Chioncel et al. (273) and Drchal et al. {21 A; 275). While this 
approach might be suitable for rather weakly correlated metals, i.e. at weak 
coupling [/, it does not capture the development of Hubbard bands. Therefore, 
it is not appropriate at intermediate and strong coupling. 
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5.4 Non-crossing approximation (NCA) 

The non-crossing approximation (NCA) is a conserving approximation corre- 
sponding to the resummation of a particular class of diagrams (those that do 
not cross). It can also be formulated as a resolvent perturbation theory in the 
hybridisation A(a;) of the Anderson impurity problem (276; 277). Hence, it is 
reliable at strong coupling where, for all orbital matrix elements, A^j(a;) <C U. 

Let us start by rewriting the non-interacting Green function of the Anderson 
model (Eq. [27]) in a form consisting of a local and a hybridisation part (bold 
symbols indicate orbital matrices): 



Then, the Dyson equation [31] of the Anderson impurity model which connects 
the non-interacting Green function G^{co) of the Anderson model and the local 
interacting Green function G{io) reads 



Since NCA is a perturbation theory in terms of the hybridisation function 
A (a;), the first step is to diagonalise the atomic problem (at) as in the Hubbard 
approximations of Section 5.2. In Section 5.2 we did so for the restriction to 
degenerate orbitals and a single Coulomb parameter U. Let us here diagonalise 
the atomic problem by formally rewriting the local Hamiltonian in terms of 
the local eigenvectors \a) and eigenenergies Eq,- 



[C?°(a;)]-^=a; + /x-e^*- A(u;) 



(101) 



where we identify the local part by 




(102) 



[G(a;)]-^ =uj + ii-e^^- A{cj) - S(u;). 



(103) 




(104) 



Ima 



Imno cc' 




(105) 



a 



In this eigenbasis, the creation and annihilation operators read 




(106) 




(107) 



a,/3 
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Here, the eigenstates denoted by |a) have one electron more than the 
states. 

A particular advantage of NCA is that this approach allows us to take into 
account the full Coulomb matrix plus the spin-orbit coupling. Of course, one 
will in practice reduce the matrix Ui^no to two-orbital terms as in Hamiltonian 
[9] of Section 4, but the spin-flip contribution of the Hund exchange coupling 
and the pair-hopping term can be included. 

The key quantity for the resolvent perturbation theory is the resolvent 
R{ijj) = {u — H)~^. Without hybridisation, this resolvent can be simply ex- 
pressed in the eigenbasis of H^^: 

^'*H = E^^I")H- (108) 

a " 



If we now take into account the hybridisation A(a;), we can define a self 

energy correction due to this hybridisation. Denoting this self energy by S{lo), 
the Dyson equation for the resolvent reads: 

Rico) = i?^*(w) + R''\u})S{lu)R{u}). (109) 



Following Keiter and Kimball (276) and Bickers et al. (277), we now express 
S{u!) in a power series in terms of the hybridisation and cut this power series 
after the lowest order term (this is the non-crossing approximation): 



mcr 



«■ J™ „./ /Of *f 



Im a'p' ' 

Here, Sa/si^^) and Ra/3{to) are the matrix elements of S{uj) and R{ijj), respec- 
tively, in the atomic eigenbasis; f{e) = 1/[1 + exp(— /3e)] denotes the Fermi 
function, and r^(e) = — 9m {A^(e + fO"*")} is the imaginary part of the 
hybridisation. The resolvent R{io) in Eq. [110] has to be determined self- 
consistently in order to yield a conserving approximation for the Anderson 
impurity model. If we expressed Eq. [110] diagrammatically we would see that 
no conduction electron lines cross. 
To employ NCA in the DMFT context, we have to determine the key quan- 
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tity of DMFT, i.e. the local Green function G{u). This is achieved via 

a,a' 13,13' 

(111) 

with the atomic partition function 

a 

With G(u;) from Eq. [Ill], we can continue the DMFT cycle, defining a new 
auxiliary Anderson impurity problem, solving it by NCA etc. 

The particular advantages of NCA are that it is a (computationally) rel- 
atively fast method — at least if not too many orbitals are involved — , and 
a conserving approximation to the Anderson impurity model. It is good at 
strong coupling and at temperatures above the effective Kondo temperature 
of the Anderson impurity model. The drawbacks are that NCA is known to 
violate Fermi liquid properties at low temperatures (below the Kondo temper- 
ature) and whenever charge excitations become dominant (278; 85). Hence, in 
some parameter ranges it fails in the most dramatic way and must therefore 
be applied with considerable care (85). The NCA has been used intensively as 
a DMFT solver, in particular, by Pruschke and coworkers. In the context of 
LDA-FDMFT, Zolfl et al. Ref. (279; 280) studied LaTiOa and cerium by NCA. 

The same shortcomings arc principally also true for the one-crossing ap- 
proximation (281; 282) which goes beyond Eq. [110] by taking into account 
additional diagrams with a single line crossing, a vertex correction. Haule et al. 
(283) used this extension recently to study the optical conductivity of cerium. 
The conserving T matrix approximation (CTMA) includes further diagrams 
beyond NCA. Kirchner et al. (284) extended this approach to calculate dy- 
namical properties which are necessary for a DMFT solver. This approach 
fulfils the Fermi liquid properties at low temperatures and might improve on 
most deficiencies of the NCA. Hence, this approximation might be suitable 
for DMFT, in particular, if the deficiency of not yielding the exact unitary 
limit, the exact height of the spectral function at a; = 0, is overcome (284). 
Most recently, also a direct expansion of the Green function in terms of the 
hybridisation has been proposed (285). 

5.5 Quantum Monte Carlo (QMC) simulations 

In the previous Sections we have already introduced several methods to cal- 
culate the Green function of the Anderson impurity model (Eq. [27]). This 
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is necessary for the DMFT self-consistently cycle (flow diagram Fig. 8), and 
quantum Monte Carlo (QMC) simulations allow to do so in a numerically 
exact way. For the Anderson impurity model, the QMC algorithm of Hirsch 
and Fye (286) is the efficient, well established approach. Hence, we will briefly 
review this QMC algorithm and a projective variant in the next two Sections. 

Let us mention however that there has been a rapid development in the 
field of QMC algorithms for the Anderson impurity model most recently. In 
particular, continuous time QMC algorithms have been developed (287; 288; 
289; 290; 291), based on a series expansions starting either from the non- 
interacting problem (287; 288) or the atomic limit (289) . Also a combination of 
series expansion and Hirsch-Fye algorithm is possible (292). One of these new 
algorithms might superseed the hirthto standard Hirsch-Fye QMC approach 
in the future. In the present state of flux, it is not clear however which of the 
new algorithms or furthergoing ones will prevail. 



5.5.1 Hirsch-Fye algorithm. In essence, the QMC technique maps the in- 
teracting Anderson impurity problem (Eq. [27]) onto a sum of non- interacting 
problems where the single particle moves in a fluctuating, time-dependent 
field. This sum is evaluated by Monte Carlo sampling, see the flow diagram 
Fig. 17 for an overview. 

Trotter discretisation 

Let us now discuss this approach in more detail. In a first step, the imaginary 
time interval [0, f3] of the functional integral Eq. [27] is discretised into A 
steps of size At = (5/A, yielding support points tx = AAr with A = 1 ... A. 
Using this Trotter discretisation, the integral Jq dr is transformed into the 
sum Ar and the exponential terms in Eq. [27] can be separated via the 

Trotter-Suzuki (294) formula for operators A and B 



A 




-AtA^-AtB 



+ 0{At) 



(113) 



A=l 



which is exact in the limit Ar 0. The single site action A of Eq. [29] can 
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Choose random auxiliary field configuration s = {s^^^/} 



Calculate the current Green function Gcur from Eq. [125] 

(Gcur)^AinA2 = i^'^^)m\i 
with M from Eq. [119] and 

tho input [G'^U' „)-'],„„ = [G-(^v)-'],„„ +S-„(^V). 



Do NWU times (warm up sweeps) 



MC-sweep (Gcur,s) 



Do NMC times (measurement sweeps) 



MC-sweep (Gcur,s) 



G = G + Gcur /NMC 



Figure 17. Flow diagram of the QMC algorithm to calculate the Green function matrix G using 
the procedure MC-sweep of Fig. 18 [closely following (35)]. 



now be written in the discrete, imaginary time as 



A-l 



amn A,A'=0 
^ / A-l 

--At Y1 U!^^'Y.'l'''r:xrmX<*x<'X- (114) 



2 

ma,m'a' A=0 



Here, the Coulomb interactions of the LDA+DMFT Hamiltonian [34], more 
specifically the inter-orbital and intra-orbital Coulomb repulsion and the Z- 
component of the Hund exchange interaction, were for convenience put into a 
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Choose M(2M — 1)A times an auxiliary field index {Xmm! era'), 

Sncw = S except for this index for which {Snew)xmm' = ~^\mm'- 



Calculate flip probability V{s Snew) = min{l, P(snew)/-P(s)} 

with P(snew)/-P(s) = det M^^— / det M° 
and M from Eq. [119]. 



res 

mn 




Figure 18. Procedure MC-sweep using the MctropoHs(293) rule to change the sign of s' 
recalculation of Gcur, i.e. tlic matrix M of Eq. [119], simplifies to 0{A^) oijciations if only one 



'^'^ ,.The 

Amm 



, changes sign. Then, Eqs. [134] and [135] can be employed (closely following (35)). 



unified form with 



^mm' 



for m = m! and a = a' 

U for m = m' and a ^ a' 

V — J for m / m' and a = a' 

V for m 7^ m' and a a' 



(115) 



Also note, that the first term was Fourier-transformed from Matsubara 
frequencies in Eq. [29] to imaginary time. 



Hubbard- Stratonovich transformation 

In a second step, the M{2M — 1) interaction terms (M denotes the number 
of interacting orbitals) of the single site action A (Eq. [114]) are decoupled by 
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introducing a classical auxiliary field s' 





;m'(CA*CA-<'A*<'A)}>(116) 



where cosh.{J""l^^,) = exp(ArC/^'^,/2). Since we need only one field for ev- 
ery pair of orbitals [{ma) and {m'a')] we can restrict tlic s^^^, sum to one 
term per pair wliich is indicated by the prime. Tliis so-called discrete Hirsch- 
Fye-Hubbard-Stratonovich transformation can be applied to the Coulomb 
repulsion as well as to the z-component of Hund's rule coupling, all in- 
cluded in C/mm'- One limitation of QMC is that it is very difficult to deal 
with terms which do not have this density-density type of form. In par- 
ticular, a Hubbard-Stratonovich decoupling of the spin-flip term of Hund's 
rule coupling leads to a 'minus-sign problem', sec Rcf. (149). Therefore, 
this spin-flip term and the pair-hopping term, i.e. the second line of the 
LDA-f-DMFT Hamiltonian [34], are usually neglected. In the particle-hole 
symmetric case, another decoupling scheme which includes the spin-flip term 
is possible without 'minus-sign problem', see Ref. (295). More recently, also 
several new algorithms (296; 297; 287; 288; 289; 290; 291) have been intro- 
duced which include the spin-flip and pair-hopping terms. This is done by new 
kinds of Hubbard-Stratonovich transformations (296; 297) or series expansions 
(287; 288; 289; 291; 292). Hence, LDA-f-DMFT (QMC) calculations with the 
correct symmetry of the interaction should be possible in the future; flrst steps 
have already been taken (291). 

By means of Eq. [116], we replace the interacting system by a sum of 
AM(2M — 1) auxiliary fields s^^j^:- This allows us to solve the functional 
integral by a simple Gauss integration since the Fcrmion operators only en- 
ter quadratically, i.e. for a given configuration s = {s^^^/} of the auxiliary 
fields the system is non-interacting. The quantum mechanical problem is then 
reduced to a matrix problem 






mAi nAg 



detM^^^. (117) 



A m'<T',m"(T" 
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Here, 

/ / 

-^ = E E E detM- (118) 

A m'a',m"a" s"^'^",;^, =±1 

is the partition function, the prime indicates that every distinct term is counted 
only once, the determinant includes the orbital, imaginary time and spin index, 
and MJ^yj is the following matrix in the imaginary time indices: 

M^^„ = Ar2[(G--^)^n + ^Ue-"^- + 1 - e-*^"^<^mn. (119) 
The elements of the imaginary time matrix are, in turn, given by 

J^mW = ""^AA' ^ '^mm'^rmn'^Xmm'J (120) 
m'cr' 

where ^^^z = 2@{a' — a + S^a' [m' — m] — 1) changes sign if {ma) and {m'a') 
are exchanged. For more details and a derivation why the matrix M enters in 
Eq. [119], see Refs. (286; 299; 56). 

Monte Carlo importance sampling 

Since the sum in Eq. [117] consists of 2^*^(2M-i) addends, a complete 
summation for large A is computationally impossible. Therefore, the Monte 
Carlo method is employed, which is often an efficient way to calculate high- 
dimensional sums and integrals. In this method, the integrand F{x) is split 
into a normalised probability distribution P and the remaining term O: 

J dxF{x) = J dxO{x)P{x) = {0)p (121) 

with 

JdxP{x) = l and P{x)>0. (122) 

In statistical physics, the Boltzmann distribution is often a good choice for 
the function P: 



P{x) = ^exp[-/3E{x)]. 



(123) 
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For the sum of Eq. [117], this probabihty distribution translates to 

P(s) = -idetM- (124) 

with the remaining term 

0(s)UnA.= [(M^^)-^]^A.nA.- (125) 

Instead of summing over all possible configurations, the Monte Carlo sim- 
ulation generates configurations xi according to the probability distribution 
P{x) and averages the observable 0{x) over these Xi. Therefore the relevant 
parts of the phase space with a large Boltzmann weight are taken into account 
to a greater extent than the ones with a small weight, coining the name im- 
portance sampling for this method. For M statistically independent addends 
Xi drawn according to the probability P{x) [this is indicated by Xi G P{x) 
in the following equation] , one gets by virtue of the central limit theorem the 
following estimate 

^0)P = Jf E 0{xi)±^^mp-{0)l. (126) 

i=l 

Here, the error and with it the number of needed addends N is nearly in- 
dependent of the dimension of the integral. The computational effort for the 
Monte Carlo method is therefore only rising polynomially with the dimension 
of the integral and not exponentially as in a normal integration. The so-called 
minus-sign problem occurs if the error [given by the variance in Eq. [126]) is 
large in comparison with the mean value. This happens particularly if 0{x) 
has contributions with positive and negative sign which almost cancel. 

In order to pick configurations x with the proper probability P{x), a Markov 
process is employed. Was x realised, a new configuration y is accepted with 
probability 

V{x^y) = mm{l,P{y)/P{x)} (127) 

Since this transition probability by Metropolis et al. (293) fulfils the detailed 
balance 

P{x)r{x ^y) = P{y)r{y ^ x), (128) 
it guarantees that the series of configurations x obey the probability P{x). 
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In our case, the probability ratio for deciding on the acceptance of a new 
auxiliary field configuration Snew is given by 

det M°"^"™ 

P(Snew)/P(s) = . (129) 

Single spin-flip updates 

A very efficient algorithm is obtained by considering a new configuration 
of the auxiliary field (spins) Snew which differs from the old one by only one 
component {Xmm' aa'): Snew = s except for (snew)Amm' = ~^'\mm'- '^^^^ 
dramatically reduces the effort to calculate the probability ratio which then 
only depends on three numbers: 

-P(Snew)/-P(s) = Rmm^-^7n/m'(7' ~ <^c^o-' -^mm'CT (1^0) 

where 

= 1 + (1 - [{M'^n-'Umx) ( e^-^-"^' ^Trl^', -!), (131) 

= 1 + (1 - [{M'^'T'Uxm'x) (e-^ - 1), (132) 
KlZ^ = mn-'Um'x{e-^ - mM'^n-'Uxmxie^ - 1). (133) 

Similarly one can show with some matrix algebra that the recalculation of 
^mn™ requires C'(2M^A^) computational operations in two steps (without 
off-diagonal orbital-elements this reduces to 0{2 x A^)). The first step is 

[{M^^-r'Uynx" = mn-'Uynx" + 

x([(M-^)-l]^;,,^A - SyxSrhm)[iM''n~XxnX"-im 

With the M^"°" from the 1. step, R^m°m'cr' recalculated according to Eq. [132] 
and the second update incorporates the change because of m': 

[(M'^'^-)-^]^A'nA" = [(Mf ^-)-l]^,.^,. + 

x([(M{^'^-)-l]^Vm'A-<5A'AW)[(M-'^-)-l]^,AnA''.(135) 

Note that these equations hold if G{t > 0) is negative, as in the standard 
definition for the Green function, and if G(0) = G{t = 0—). Often in QMC, 
the Green function is defined differently, i.e. with the opposite (positive) sign. 
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Since roughly M(2M — 1)A single spin flips have to be tried before arriving 
at a truly independent Snew) the overall cost of the algorithm is 

2M(2M - l)M'^h? X number of MC-sweeps (136) 

in leading order of A. It reduces considerably to 

2M(2M - 1)A^ X number of MC-sweeps (137) 

if orbital off-diagonal elements are zero. 

The advantage of the QMC method (for the algorithm see the flow diagrams 
Figs. 17 and 18) is that it is (numerically) exact. It allows one to calculate 
the one-particle Green function as well as two-particle (or higher) Green func- 
tions. On present workstations the QMC approach is able to deal with up 
to seven interacting orbitals and room temperature or higher temperatures, 
for typical values of the Coulomb interaction U and the LDA bandwidth W . 
Since the QMC approach calculates G{t) or G{iuJn) with a statistical error, 
it also requires an analytical continuation to obtain the Green function G{uj) 
at real (physical) frequencies u) or the physically relevant spectral function 
A{ijj) = — ^ImG(u;), see Section 5.5.4. Very low temperatures are not accessi- 
ble in QMC because the numerical effort grows as A^ ~ 1/T^. For these low 
temperatures or zero temperature, a projective QMC method was developed 
recently. We will discuss this variant in the next Section. 



5.5.2 Projective quantum Monte Carlo (PQMC) simulations for T — 0. 
Often interesting many-body physics occurs at low temperatures. In this case, 
the Hirsch-Fye (286) QMC algorithm, which was introduced in the previous 
Section, is not applicable as discussed above. For lattice QMC simulations an 
alternative projective quantum Monte Carlo (PQMC) method was developed 
by White et al. (300). This PQMC algorithm converges, according to Assaad 
and Imada (301), faster to the groundstate than the finite temperature lattice 
QMC algorithm of Blankenbecler et al. (302). The general idea of the PQMC 
is to start with a trial wave function 1^*0) and to project onto the ground state 
I^Gs) via 



I^Gs) = Hm '/^ :■ (138) 

Y'(^'o|e-^-^|*o) 



0-^00 



If l^'o) has any overlap with the (unique) ground state only this state with 
the lowest energy and, hence, the largest contribution to e~^^ will prevail in 
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the limit ^? ^ oo. In the DMFT context, one often calculates the Anderson 
impurity model in the thermodynamic limit with an infinite number of bath 
sites A?^ — oo. In this case, ground state and low-lying excited states have 
the same energy to leading order in and yield the same Green function. 
Hence, we can project onto such low-lying excited states and obtain the correct 
Green function and other expectation values, so that Anderson's orthogonality 
catastrophe (303; 304; 305) is irrelevant for the PQMG algorithm, as was shown 
by Feldbacher et al. (306). 

With the ground state given by Eq. [138], zero temperature observables can 
be calculated as follows: 

{O) = {^GS\ O |*GS) = lim ' ' ^ , - ,^ ^ °^ (139) 

For a finite value of 9, this expectation value can be calculated by (projec- 
tive) QMC simulations which decouple the projector e~^/'^^ by Hubbard- 
Stratonovich transformations, after discretising the imaginary time. 

Let us now focus on the Anderson impurity model for which Feldbacher et 
al. (189) developed a new PQMC method. This PQMG is related to the one of 
White et al. (300) for lattice many-body systems in a similar way as the QMC 
of Hirsch and Fye (286) to that of Blankenbeclerei al. (302) for lattice QMG. 
The biggest difference is that, for the Anderson impurity model, one is dealing 
with matrices in the imaginary time index, instead of matrices in the lattice 
indices for the lattice QMC algorithm of Blankenbeclerei al. (302): The Hirsch- 
Fye algorithm is directly formulated in terms of the non-interacting Green 
function of the Anderson impurity model . Within the DMFT iteration 
scheme, we usually do not even define a lattice for the auxiliary Anderson 
impurity model anymore. This leads to the difficulty of how to define the trial 
wave function which for a lattice problem is a more straightforward task. 

Feldbacher et al. (189) overcame this difficulty by considering instead of Eq. 
[139] an artificial finite temperature problem: 

= ^ gfy OH • ^1^°) 

p—*00 'JYg— /5^i0 g — P-H AIM 

In the limit /3 — oo, this projects onto the ground state of the Hamilton 
operator ^o- If this ground state is |*o) and has energy Eq we have 

lim Tre-'^^°d e-^^°(*o| O |*o)- (141) 

^— »oo 



Hence, the finite temperature problem of Eq. [140] becomes equivalent to the 
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projection of Eq. [139]. 
As a special case of Eq. [140], the ^-projected Green function is given by 

G,„.(T,r') = -(rq(r)4(r')) (142) 

(^o|re-^^-"Q(T)cL(r')|^o) 



(143) 



Here, T is the Wick time ordering operator for the r's and 0{t) = e^^Oe~'^^ 
as before. 

A natural choice for the trial wave function is the ground state of a non- 
interacting Anderson impurity model [[/ = in Hamiltonian [26]). As was 
shown by Feldbacher et al. (189), one then arrives at the same algorithm as 
the Hirsch-Fye finite temperature algorithm. The only difference is that in Eq. 
[119], instead of the finite temperature {Q^)~^ = + S, a zero temperature 
{Q^)~^ enters as a A x A matrix. While at half- filling the finite temperature 
^^j„(t, 0) decreases (in terms of the absolute value) from T = 0toT = /3/2 
and increases again from r = /3/2 to r = the zero temperature ^[^^(r, 0) 
continues to decrease from t = 9 (2 to t = 9. 

With these rather small changes to the QMC code, PQMC simulations are 
possible. In contrast to the finite temperature QMC, there are no thermal 
fluctuations. This is the reason why these PQMC calculations converge much 
faster to the ground state properties (as a function of 9) than the Hirsch-Fye 
QMC (as a function of (3 which corresponds to the same numerical effort). Let 
us mention that similar projective QMC versions of the continuous time QMC 
algorithms (287; 288; 289) are possible (298). 



5.5.3 Fourier transformation from r to iuj,,. In the previous two Sections, 
we have introduced the Hirsch-Fye QMC and projective QMC simulations for 
the Anderson impurity model. In the context of DMFT, this Anderson im- 
purity model is determined self-consistently, see flow diagram Fig. 8. While 
the Dyson equation [32] of this self-consistency is formulated in terms of Mat- 
subara frequencies, the QMC calculations of the Anderson impurity model 
use a discrete set of A imaginary time support points. Hence, we have to 
overcome a final technical obstacle: the Fourier transformation from r to iui, 
and vice versa. This is in particular problematic since wc have only a limited 
number A of r points, which only allow us to calculate the Green function 
for an equal number of Matsubara frequencies by means of a discrete Fourier 
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transformation. On the other hand, the tail at high Matsubara frequencies is 
important since it is responsible for the jump of the Green function at r = 0: 

Glm{T = 0+) = -Sim + Gimir = O'). 

Therefore, a discrete Fourier transformation docs not work. It would yield a 
Green function which oscillates considerably around the correct G{t). Different 
paths are used in the literature to overcome this obstacle: 

Georges et al. (56) do a spline interpolation of G(r), resulting in arbitrar- 
ily many support points and, hence, enough Matsubara frequencies. Jarrell 
(299) extend the number of Matsubara frequencies by employing the iterated 
perturbation theory result at high frequencies. 

Ulmke et al. (104) use a smoothing procedure which replaces Qmmii'^i^) by 

~0 ^'T 

This is Fourier tranformcd to imaginary time, and this Fourier transform is 
in turn used as the non-interacting Green function of the Anderson impurity 
model. After the QMC simulation yielded the output Gmm{Tx), the process 
is reversed: From the Fourier transform of GmmiTx), i.e. Gmmii^u)-, the in- 
verse of Eq. [145] yields Gmm{i^v)- The new self energy is then Emm(i<^t/) = 
(iuJu) ■ This approach generates smooth Green functions 
G{t\) and reproduces the correct At limit. 

In their implementation, McMahan et al. (229) use a constrained fit to the 
output QMC impurity Green function G{tx): 

G{T) = J2wifi{T). (146) 



The basis functions are /i(r) = — e~^*'^/(e~^^*-|-l) and have Fourier transforms 

fi{iio) = l/{iLO — £i), hence allowing for determining G{iuJu) at > A Mat- 
subara frequencies. The difference to the spline-fit is that in Eq. [146] every 
fit coefficient is determined by the local behaviour in a frequency interval, not 
by the local behaviour in an imaginary time interval. 

The fit of Eq. [146] has some additional constraints: Wi > 0, ^(0"^) is pre- 
cisely the QMC value, G(0+)+G(/3-) = -1, and ^G{Q+)+-^G{p-) = g2, where 
gm is the (za;)""* high-frequency moment of G{iu>). For the last constraint, the 
second moment g2 is obtained from the relation G~^{iLL>) = Q~^{iLL!) — T,(iu;) 
which implies g2 = 52 + sq, where these are the indicated moments of G{iLo), 
Q{iuj), and S(iLj), respectively. Note that 52 is known since Q is input to the 
QMC, and for sq the analytical high frequency behaviour (known from the 
Hartree Fock) is taken: sq = S(aj = oo) = 5]Hartrcc-Fock_ M degenerate or- 
bitals and without Hund's exchange coupling this is: sq = 2M— 1/2M[1-K?(t = 
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O"*")]?/. Using these constraints and A/4 equally spaced Cj, the agreement be- 
tween Eq. [146] with the QMC data for G(r) is optimised. 

For a faster convergence, McMahan et al. (229) also separated the self energy 
T,^^^{iu;) = AS(za;) + into a constant Hartree-Fock contribu- 

tion and a frequency-dependent rest AT,(iuj). After every QMC calculation, 
the authors did computationally inexpensive iterations during which AT,{iLu) 
was kept fix, but s^^*''^^"'^"'^'^ and the number of interacting electrons were 
adjusted self-consistently until convergence. 

For the PQMC simulations discussed in the previous Section, the Fourier 
transformation is even more challenging since the zero temperature Green 
functions extend from to oo, instead of to /? at finite temperatures T = 1/(3. 
With the Green function known only at A discrete r points, one has not 
only to interpolate between the r points as for finite temperatures but also 
to extrapolate to r = oo. To this end, Feldbacher et al. (189) employed the 
maximum entropy method (sec Section 5.5.4) which allows us to calculate from 
A support points the zero temperature spectral function A{u) at real valued 
frequencies. From this spectra function, the Green function at any Matsubara 
frequency Ui, can be easily calculated as 



5.5.4 Maximum entropy method. Since, in QMC, the Green function is cal- 
culated on the imaginary (Matsubara) axis, we have to do an analytical contin- 
uation to real frequencies for getting the spectral function A{lx)) = — -ImG(a-') 
which is of direct physical interest since it can be measured, e.g. by photoemis- 
sion experiments. Because of the statistical QMC error, the standard approach 
for doing the analytical continuation is the maximum entropy method. In the 
following the basic concept will be briefly discussed, for a detailed review by 
Jarrell and Gubernatis see (188). 

Starting point is the Fourier-transform of the spectral representation of the 
Green function 



This equation already shows that the analytical continuation is an ill- 
conditioned problem: The kernel of Eq. [148] is very small for large frequencies 
u so that large changes in A[lo) have only a small impact on G{t). 
In this problematic situation, the maximum entropy method is used. It starts 




(147) 




(148) 



— oo 
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with the entropy 

uio/2 

S = - J dw A{uj) In [A{uj) loo] (149) 

-LOo/2 

as the a priori probabihty P{A) ~ exp(TS') for a given spectrum A(Ld), with 
adjustable parameter T. Besides this constant a priori probability, A{u!) has 
to yield the QMC-calculatcd G{t). How good this is achieved can be measured 
by the usual value for the quadratic difference between given (calculated) 
G{t) and the one obtained from A{uj) via Eq. [148]. This value gives the 
conditional probability for G(r) with a given A{uj): 

P{G\A) = e---^\ (150) 

Following the Bayes theorem a priori and conditional probability yield the a 
posteriori probability for having a spectrum A{lo) if G{t) was calculated in 
QMC: 

P{A\G) ~ P{G\A)P{A) ~ e^^-2X^ (151) 

The result of the maximum entropy method is the most likely spectrum A(uj), 
maximising P{A\G). This is a well defined statistical method for doing the an- 
alytical continuation. Because of the statistical error in the QMC and because 
Eq. [148] is ill conditioned for large \uj\ the maximum entropy can however not 
resolve fine details at large frequencies, i.e. in the Hubbard band. In contrast, 
features at small frequencies, such as height and width of the central peak, 
and the overall weight and position of the Hubbard bands are reliable. 

5.6 Comparing different DMFT solvers for Lui-xSrxTiOs 

After the discussion of different methods X to solve the DMFT self-consistency 
equations in the previous Sections, let us now compare how this solver reflects 
in the LDA-I-DMFT(X) results for a specific material: Lai-a^Sra^TiOs. Such 
a comparison has been carried out by Nekrasov et al. (159), and we will 
recapitulate their findings here. The stoichiometric compound LaTiOs is a 
cubic perovskite with a small orthorhombic distortion (the distorted angle 
is Z Ti-O-Ti 155°; see Ref. (307)), and is an antiferromagnetic insulator 
below Tjv = 125 K (308; 309). Above Tjv, or at low Sr-doping x, LaTiOs is a 
strongly correlated, but — except for the distortion — simple paramagnet with 
one 3d electron on the trivalent Ti sites. Since the aim of this Section is the 
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Figure 19. Densities of states of LaTiOs calculated with LDA-LMTO. Upper figure: total DOS; 
lower figure: partial t2g (solid lines) and eg (dashed lines) DOS (reproduced from Ref. (159)). 



comparison of DMFT solvers we neglect the small orthorhombic distortion, 
i.e. consider a cubic structure with the same volume. 



LDA calculations for the cubic crystal structure 

Fig. 19 shows the LDA DOS for undoped LaTiOs. Thereby, Anisimov et 
al. (31) and Nekrasov et al. (159) approximated the crystal structure by a 
cubic one with the same volume. LDA+DMFT calculations including the 
orthorhombic distortion have been recently carried out by Pavarini et al. 
(310; 311) and Craco et al. (312). Here however, our main goal is the com- 
parison of different DMFT solvers in a realistic material calculation. For this 
intercomparison the approximate cubic structure is sufficient. 

In Fig. 19, the oxygen bands range from —8.2 eV to —4.0 eV and arc com- 
pletely filled so that Ti is three- valent. The cubic crystal field splits the Ti 
2d bands into two empty e^-bands and three degenerate t2g-bands. well sepa- 
rated from the other bands. For the low-energy behaviour it is hence possible 
to restrict ourselves to these degenerate t2g orbitals within the approximation 
introduced in Section 4.1.5, i.e. using the LDA DOS (Eq. [69]) instead of the 
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Figure 20. Spectrum of Lao.94Sro.o6Ti03 as calculated by LDA+DMFT(X) at T = 0.1 eV 
1000 K) and [7 = 4 eV employing the approximations X=IPT, NCA and numerically exact 
QMC. Inset left: Behaviour at the Fermi level including the LDA DOS. Inset right: X=IPT and 
NCA spectra at T = 80 K (reproduced from Ref. (159)). 



full one-particle Hamiltonian H^j)js^ of (Eq. [32]). We take Sr-doping x into 
account by adjusting the LDA+DMFT chemical potential to n = 1 — x = 0.94 

t2g electrons. 
Method matters 

In Fig. 20, we present the LDA+DMFT(X) spectrum of Lao.94Sro.o6Ti03, 

calculated for the impurity solvers X=IPT, NCA and QMC. Qualitatively, all 
three methods X yield the characteristic three peak structure, consisting of 
lower Hubbard band, quasiparticle peak and upper Hubbard band. 

Quantitatively, however, we see pronounced differences: The IPT quasiparti- 
cle peak is very narrow and hence not visible at high temperatures (in the main 
panel), also the shape of the IPT Hubbard bands is different. NCA is much 
closer to the (numerically) exact QMC result than IPT, but nonetheless un- 
derestimates the width of the quasiparticle peak by a factor of two. NCA also 
violates the Luttinger pinning of the spectral function (66), a known deficit 
(313) seen clearly in the left inset, and puts the lower Hubbard band to close 
to the Fermi energy. The comparison of Fig. 20 hence shows that, at least 
on a quantitative level, it matters which method is employed as an impurity 
solver for the DMFT equations. 
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6 Realistic material calculations with DMFT 

In this Section, we will review those electronic structure calculations which 
have been accomplished hitherto using DMFT for actual material calculations. 
The results are subdivided into different material classes: / electron systems 
in Section 6.1, transition metals in Section 6.2, their oxides in Section 6.3 and 
other materials in Section 6.4. Some highlights arc discussed in more detail. 
Not discussed are LDA+U calculations, which would require a review on its 
own, and model DMFT calculations without material-specific LDA or GW 
input. 

6.1 / electron systems 

The state-of-the-art LDA-I-DMFT calculations for /-electron systems take into 
account all spdf valence orbitals and all hybridizations between them, but re- 
strict the DMFT Coulomb interaction to the / orbitals in the LDA-I-DMFT 
Hamiltonian (Eq. [34]). This is a very reasonable starting point since the spd 
orbitals are much more extended and hence less strongly interacting. Whether 
the interaction of the d electrons leads to corrections, which cannot be com- 
pletely ignored, remains however an open question since presently a DMFT 
calculation with all 24 interacting d and / orbitals would be too involved, at 
least when using a more rigorous impurity solver. 

Two of the very early successes of LDA-I-DMFT were the calculations of 
the Mott transition in plutonium by Savrasov, Kotliar and Abrahams (165) 
and of the cerium volume collapse transition by Zolfl et al. (280) and Held, 
McMahan and Scalettar (154). Savrasov et al. (165; 150; 314) studied the 
5-phase of Pu which is not well described by LDA, as is evident from the 
underestimation of the Pu volume by 30%. Within LDA-I-DMFT, electronic 
correlations drive the system towards a Mott transition. With these strong 
electronic correlations and a Coulomb interaction of ?7 ~ 4eV, the volume 
of 5-Pu comes out correctly (165). Subsequently, Dai et al. (315) extended 
these calculations by including lattice dynamics, allowing for the determina- 
tion of the first LDA-I-DMFT phonon spectrum. This theoretical prediction 
well agrees with the later experiments (316). Concerning the magnetic prop- 
erties of 6 Pu, Pourovskii et al. (317) predicted the absence of dynamical and 
static magnetic moments in agreement with experiment. A detailed review of 
these activities, which are of possible relevance for nuclear waste (165), can 
already be found in Ref. (40). Hence, we will not discuss these calculations 
in more detail here. But let us mention that, more recently, also americium 
(318), PuCoGas (319), different actinide monochalcogenides (320), Pu and 
Am compounds (321), and the rare earth elements Nd and Pr (322) have been 
investigated. In Ref. (322), also the effect of the spin-orbit coupling, which 
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is certainly important for / electron systems, was taken into account. Also, 
one of the very first LDA+DMFT(Hubbard-I) calculations was for the mixed 
valence 4f material TmSe (32) . 



6.1.1 Volume collapse transition in cerium. Let us now discuss in de- 
tail elemental cerium, the material hitherto most intensively studied with 
LDA+DMFT by various groups: Zolfl et al. (280), McMahan et al. (154; 229; 
322), Haule et al. (283), and Amadon et al. (323). Under pressure, Ce un- 
dergoes the volume collapse, or a-^ transition, with a volume change of 15% 
at room temperature. This transition is isostructural within the face-centred- 
cubic (fee) lattice structure and fades away for temperatures above the critical 
point at T = 600±50K (324), for reviews see (324; 325; 326; 327). The volume 
collapse reflects in a dramatic transfer of spectral weight: In the a phase there 
is a large peak at the Fermi energy; whereas this peak is very much reduced in 
the 7 phase, albeit there is still some spectral weight at the Fermi energy as 
to be expected for a metal. In accord with these findings, the optical conduc- 
tivity is higher in the a phase where the frequency dependent scattering rate 
is characteristic for a Fermi liquid behaviour with an effective mass of about 
20TOe, see Ref. (328). Concerning the magnetic susceptibility, a Ce behaves 
like a Pauli paramagnet at room temperature but has a Curie- Weiss form at 
higher temperatures and in the 7 phase (324). Despite these dramatic differ- 
ences in the spectrum, optical conductivity, and susceptibility, the number of 
4/ electrons does not change significantly and is close to one across the 01-7 
phase transition, as was revealed by Myon decay experiments (324). 

Although some alternative theories have been proposed, see e.g. (329; 330; 
331), the general belief is that the a-7 transition has an electronic origin. The 
first electronic theory, the promotional model (332; 333), which assumed the 
electronic configuration to change from 4/°(spd)^ in a Ce to 4f^{spd)^ for the 
7 phase, was however dismissed since it is at odds with the above mentioned 
Myon experiments and also with LDA calculations (334). Instead, Johansson 
(334) proposed a Mott transition (MT) model for the 4/ electrons, assum- 
ing the 4/ electrons in the a phase to be itinerant whereas they are localised 
(Mott insulating) in the 7 phase which has a reduced LDA 4/ bandwidth. 
Of course, 7 Ce as a material remains metallic due to the other {spd) elec- 
trons. In subsequent efforts to treat this scenario within LDA, Johansson (335) 
employed standard LDA calculations for the a phase and treated the 4/ elec- 
trons in the 7 phase as localised spins. Qualitatively, this yields similar results 
as the more sophisticated self-interaction corrected LDA calculation and the 
LDA+U method which were later performed (336; 337; 338; 339; 340) along 
with orbitally polarised calculations (341; 342; 336). 

At the time Johansson was working out the MT scenario, another phe- 
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nomenon of electronic correlations was finally understood: the Kondo effect. 
Based on the physics of the Anderson impurity model, Allen and Martin 
(343; 344) proposed a Kondo volume collapse (KVC) scenario for the a-7 tran- 
sition which appeared to be in conflict with the MT scenario. Both pictures 
agree that, at the experimental temperatures, the large volume 7 phase has 
strongly correlated (localised) 4/ electrons with a 4/^ moment and a Curie- 
Weiss magnetic susceptibility. But they differ considerably for a Ce: The MT 
scenario of Johansson assumes a weakly correlated (itinerant) a phase, with a 
'single' peak at the Fermi energy in the 4/ spectrum as on the left hand side 
of Fig. 2, only with some additional structures. As in the LDA, this itinerant 
phase has no 4/^ moment at all. The KVC scenario on the other hand assumes 
continued strong correlation in the a phase with a three peak structure, in- 
cluding a central Abrikosov-Suhl resonance and two side peaks at considerably 
higher and lower energies (345), as in the middle of Fig. 2. While the MT 
scenario envisages the coexistence of two phases and a first order transition 
between these, the KVC assumes the /-valence hybridisation to increase upon 
pressure so that the Kondo temperature which depends exponentially on this 
hybridisation changes dramatically. This leads to a correlation contribution to 
the energy (which is roughly proportional to the Kondo temperature) with a 
negative curvature as a function of volume. Therefore, one can lower the energy 
with a mixed low and high volume phase by a Maxwell (tangent) construction 
(343; 344), indicating a first order transition similar to the vapour-liquid tran- 
sition. Attempts to do realistic material calculations on the basis of this KVC 
scenario started with Ref. (345) and were later continued (346; 347). In these 
calculations, the /-valence hybridisation is determined from LDA followed by 
a many-body treatment of the Anderson impurity model with the seven 4/ 
orbitals. 

Studying the Hubbard model and the periodic Anderson model^ which are 
the simplest models for the MT and KVC scenario, respectively. Held et al. 
(348; 349) more recently stressed the similarity of these two models within 
DMFT and, hence, the MT and KVC scenarios. Actually, one can integrate 
out the valence electrons of the periodic Anderson model, resulting in an ef- 
fective one-orbital model such as the Hubbard model but with a frequency 
dependent kinetic energy term. Hence, it might not be surprising that also the 
physics of the two models is very much the same: The spectrum shows a three 
peak structure with a lower and an upper Hubbard band and a central quasi- 
particle resonance in between. The local moment of the a phase is screened 
at low energies for both the Hubbard and the periodic Anderson model. This 
shows that the differences between the MT and the KVC scenario are due to 



^Whilc the KVC is based on the Anderson impurity model, a more reaUstic treatment should start 
with its periodically extended version since the 4/ electrons on every Ce site interact. 
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the LDA treatment of the a phase: In the a phase, the MT scenario misses 
electronic correlations and, therefore, the three peak structure of the spec- 
trum. While many physical aspects of these models are similar, there are also 
some noteable differences. In particular, there is no Mott-Hubbard transition 
in the periodic Anderson model at zero temperature, even though the finite 
temperature behavior is very similar to that of the Hubbard model including 
a region of two coexisting phases, i.e. metal and insulator (350). 

Following these model studies, realistic LDA+DMFT calculations provided 
for an accurate ab-initio description of Ce. Fig. 21 (left panel) shows the evo- 
lution of the k-integrated 4/ spectral function with increasing volume. At a 
very small volume, = 20 A'^, most of the spectral weight is in a central quasi- 
particle peak or Abrikosov-Suhl resonance at the Fermi energy, similar as in 
the LDA. But, Hubbard side structures are already discernible, indicating that 
there are already electronic correlations, albeit these are not yet extraordinar- 
ily strong. With increasing volume, more and more spectral weight is trans- 
ferred to the Hubbard side bands; electronic correlations increase. Approach- 
ing the experimental volumes of the a-7 transition, which occurs between 28.2 
and 34.4 A^ at room temperature, the three peak structure becomes much 
more pronounced: We see a sharp quasiparticle resonance at F = 29 A^. From 
the experimental a to the experimental 7 phase volume (from V = 29 to 
34 A^ in Fig. 21), the weight of the quasiparticle peak shrinks dramatically and 
fades away at even larger volumes. This large volume phase with a two peak 
structure can also be described by LDA-l-DMFT(Hubbard-I) and LDA-I-U 
calculations (not shown). Across the transition from 7 to a Ce, also the local 
magnetic moment (m|) is reduced by 5% (229) but not lost, in contrast to 
the Mott transition scenario of Johansson (335). 

Adding the spd valence spectrum to the 4/ spectrum and multiplying with 
the Fermi (inverse Fermi) function, the total LDA-I-DMFT spectrum is com- 
pared to the photocmission spectroscopy (PES) spectrum (351) below the 
Fermi energy and to the Bremsstrahlung isochromatic spectroscopy (BIS) 
(352) spectrum above the Fermi energy in Fig. 21 (right panel). The agree- 
ment between theory and experiment is very good. Note, that there are no 
free parameters in the LDA+DMFT(QMC) results of McMahan, Held and 
Scalettar (353; 229) since the /-electron Coulomb interaction and the double 
counting correction have been determined by constrained LDA calculations 
and the experimental resolution has been taken from (280). Particularly good 
is the agreement of the quasiparticle peak around the Fermi energy for both 
a and 7 Ce, but also the position of the upper and lower Hubbard bands 
are approximately correct. The biggest differences can be found in the upper 
band which is broader in experiment and has some inner structure. It was 
argued (280; 229) that this is due to the Hund exchange interaction which was 
neglected in both, LDA-FDMFT(NCA) (280) and LDA-FDMFT(QMC)(229), 
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Figure 21. Left: Evolution of the 4/ spectral function A{u;) with volume at T = 632 K. The curves 
at different volumes axe shifted as indicated by the base line. When going from small to large 
volumes, the weight of the central quasiparticle peaJc is dramatically reduced at volumes 
corresponding to the experimental a-7 transition from V = 29 to 34A'^. Right: Comparison of the 
parameter- free LDA-|-DMFT(QMC) spectrum with experiment (circles), as collected in Ref. (345) 

(reproduced from (229)). 



calculations. This exchange interaction has only a minor effect for the occu- 
pied states since these consist mainly of single occupied configurations. But, 
it is important for the upper Hubbard band because the Hund exchange splits 
these doubly occupied states into multiplets. 

With decreasing volume, we saw the development of a quasiparticle peak 
in the 4/ electron spectrum. This quasiparticle physics is associated with an 
energy gain which cannot be captured by LDA and which gives rise to a cor- 
relation energy with a negative curvature (353; 229) at low temperatures (not 
shown). Hence, it also reflects in the total energy shown in Fig. 22: At high 
temperatures and for the LDA -l-DMFT (polarised HF) results, we see a simple 
-E'tot(^) curve with a single minimum, the equilibrium volume. But at lower 
temperatures, the negative curvature of the correlation energy leads to a side 
structure. Since all energy contributions except for the correlation energy have 
a positive curvature, the negative curvature of the correlation energy is largely 
compensated, and we see a very shallow region at T = 0.054 eV. Within the 
numerical error bars, it is difficult to decide whether we already have a nega- 
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Figure 22. Left: Total LDA+DMFT{QMC) and LDA-|-DMFT(polarised HF) energy as a function 

of volume at three temperatures. The negative curvature of the correlation energy results in the 
development of a side structure, visible as a deviation from the LDA+DMFT(polarised HF) energy. 
The long dashed line is the curve which corresponds to the pressure of the a-7 transition: 
E = — PexpV (reproduced from (229)). Right: Results for the total energy from (323), including 
also the experimental free energy, (reproduced from (323)). 



tive curvature or whether shghtly lower temperatures are needed. A negative 
curvature of the total energy will give rise to a Maxwell construction and a first 
order phase transition. The region in which the LDA-I-DMFT(QMC) results 
show these tendencies are in agreement with the experimental a-7 transition 
which is marked by the arrows in Fig. 22. The slope of the shallow minimum is 
also consistent with an experimental pressure of —0.6 GPa (long-dashed line) 
given by the a-7 transition pressure extrapolated to T = (335). 

Since vertex corrections do not contribute to the optical conductivity within 
DMFT (85; 86), the optical conductivity can be calculated directly from the 
DMFT G(k, Ld) and the dipole matrix elements. Haule et al. (283) calculated 
these dipole transition matrix from the LDA wave functions and used the 
one-crossing approximation (OCA) for solving the auxiliary DMFT impurity 
model. The resulting optical conductivity in Fig. 23 shows shows a clear Drude 
peak for the a phase, while for the 7-phase, the conductivity is much smaller 
for uj ^ and no Drude peak is discernible. The basic features of the the- 
oretical optical conductivity agree with experiment (lower panel of Fig. 23). 
Hence, we can altogether conclude that LDA+DMFT correctly describes the 
thermodynamic and spectral properties of the a-7 transition, as well as the 
4/ occupation of n/ ~ 1 in the vicinity of the transition (229). 
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Figure 23. Comparison of LDA+DMFT optical conductivity (top panel; Ref. (283)) with 
experiment (lower panel; Ref. (328)). The LDA+DMFT spectra are for 580 K (a-Ce) and 1160 K 
(7-Ce); the experimental temperature was 5 K (a-Ce) and 300 K (7-Ce), respectively, (reproduced 

from Ref. (283)). 



6.2 Transition metals 

6.2.1 Ferromagnetism in Fe and Ni. Electronic structure calculations for 
transition metals and transition-metal-transition-metal alloys so-far concen- 
trated mainly on the ferromagnets Fe and Ni, starting with early LDA-I-DMFT 
calculations by Drchal et al. (274; 275) and Lichtenstein et al. (158). Since 
some aspects of ferromagnetism in Fe and Ni, in particular the ferromagnetic 
moment, are well described by conventional LDA calculations, the question is: 
Are transition metals strongly correlated as the importance of the 2>d orbitals 
suggest or not? In other words: Is an LDA-I-DMFT calculation necessary for 
iron and nickel or is LDA sufficient? 

Certainly not described by LDA is the famous -6eV satellite in Ni. Using 
LDA-I-DMFT, Lichtenstein, Katsnelson and Kotliar (158; 36) reinvestigated 
this element, and did indeed find a satellite at about -6eV, see Fig. 24 (left 
panel). This spectral feature could hence be explained by LDA-I-DMFT as a 
Hubbard band in the majority-spin spectrum. Later, this finding was also con- 
firmed by GW-I-DMFT calculations by Biermann, Aryasetiawan and Georges 
(42; 43; 44), reporting similar results as LDA-I-DMFT for Ni, see Fig. 24. 

Along with the satellite peak, Lichtenstein et al. (158; 36) found for the 
paramagnetic phase a Curie susceptibility which indicates the presence of lo- 
cal (unordered) magnetic moments of size 3.09 and 1.50 for Fe and Ni, 
respectively. This Curie susceptibility is in agreement with experiment, see 
Fig. 25. Along with the satellite, the local moment is clearly a correlation 
effect, since within L(S)DA the local magnetic moment fades away with the 
magnetisation Similarly, LDA-I-DMFT calculations for a and 7 Ce yield lo- 
cal magnetic moments in the paramagnetic phase. The absolute value for the 



February 2, 2008 18:28 Advances in Physics ESCwDMFT'16 



Electronic structure calculations using dynamical mean field theory 



85 




-4 -2 

Energy (eV) 



5 
Energy (eV) 



Figure 24. Spectrum (k-integrated) of Ni [left; LDA+DMFT (solid lines), spinpolarised LDA 

(LSDA) (dotted lines); right: GW+DMFT]. The two lines represent the minority- and 
majority-spin spectrum respectively. At roughly -6eV, a satellite peak is clearly visible in the 
majority-spin spectrum (reproduced from (158) and (42), respectively). 




Figure 25. Magnetisation below and magnetic susceptibility above the Curie temperature for iron 
(squares) and Ni (circles). The LDA-I-DMFT results (open symbols) are compared to experiment 
(full symbols) (reproduced from Ref. (158); experiments from Ref. (354; 355)). 



Curie temperature is somewhat overestimated, to a lesser extent in Ni (10%) 
than in Fe (80%). This is (i) because the DMFT neglects non-local correlations 
such as spin waves and (ii) because the LDA+DMFT calculations employed 
a Z2-symmetric Hund's exchange instead of a SU(2)-symmetric so that trans- 
verse spin fluctuations are suppressed. Let us add that the general features of 
the susceptibility in Fig. 25 have been found for the ferromagnetic transition 
of a simple model: the single band Hubbard model on a fee lattice (100). 

More recently, Minar et al. (356; 357) employed the Korringa-Kohn-Rostoker 
(KKR) method for LDA+DMFT studies of Fe, Co, Ni and the alloy Fe^^Nii.^;, 
describing -among others- the Fano effect in these ferromagnets, and Grechnev 
et al. (358) studied Fe, Ni and Co surface spectra. Ferromagnetic Ni was also 



February 2, 2008 



18:28 Advances in Physics ESCwDMFT'16 



86 K. Held 

at the centre of the reahstic Gutzwiller calculations by Biinnemann et al. (218; 
219; 220); for a summary of these results see Ref. (221) and, supplementarily, 
Ref. (359) for the angular resolved spectrum and Ref. (360) for the total energy. 
LDA+DMFT calculations for magnetic multilayers of transition metals have 
been carried out recently by Chioncel et al. (273; 361); also note Section 6.4.1 
discussing half-metallic ferromagnetism in Heussler alloys. 

6.3 Transition metal oxides 

Transition metal oxides are an ideal laboratory for the study of electronic cor- 
relations in solids, showing a rich spectrum of physical phenomena, ranging 
from the Mott-Hubbard transition in V2O3, to high temperature supercon- 
ductivity (among others in Sr2Ru04) and the heavy Fermion behaviour in 
LiV204. The materials mentioned and others have been studied meanwhile 
by LDA+DMFT. In these materials, the 3(i bands are comparatively narrow 
with width W 2 — 3 eV so that electronic correlations, induced by the 
local Coulomb interaction U k, 3 — 5 eV, are strong. Hence, in Fig. 2, tran- 
sition metal oxides are neither in the weakly correlated region {U /W <C 1) 
nor is the opposite limit {U /W ^ 1) appropriate. These materials are in the 
'in-between' regime, U /W = 0{1). If the ?>d transition metal oxide is metallic 
this phase is strongly correlated with a quasipartide peak at the Fermi energy 
and Hubbard side bands. But depending on the material, also Mott insulating 
behaviour and the Mott-Hubbard metal-insulator transition can occur. 

In all LDA+DMFT calculations for transition metal oxides, the 
LDA+DMFT Hamiltonian (Eq. [34] ) was restricted to the low-energy orbitals 
at the Fermi energy, typically the t2g or eg orbitals. These Wannier orbitals 
represent mixtures between mainly the transition metal d and the oxygen 
p orbitals. So far, LDA+DMFT calculations have not succeeded in taking 
into account a larger basis including oxygen p orbitals. The naive inclusion 
of non-interacting oxygen orbitals result in strong deviation from an integer 
occupation of the interacting d orbitals. Consequently, electronic correlations 
are too weak (362). This should be overcome if the p-d and the p-p interactions 
are included. However, such calculations arc computationally very demanding 
presently. Often, also the simplification for transition metal oxides which we 
discussed in Section 4.1.5 and which allows us to do the DMFT calculation 
with the LDA DOS only was employed. For non-cubic systems this is an ap- 
proximation, which is however very reasonable as long as the crystal is not too 
strongly distorted from the cubic symmetry. Typically, every transition metal 
ion is still surrounded by an octahedron of oxygen ions. However, this octa- 
hedron is then not perfect anymore, but tilted and distorted. Consequently, 
there are orbital off-diagonal elements between the Cg and t2g orbitals in the 
LDA. 
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In the following, LDA+DMFT calculations for some transition metal oxides 
will be presented in detail. Besides these studies, LDA+DMFT has been also 
applied to Cr203 for which Craco, Laad and Miiller-Hartmann (475; 476) 
analysed orbital correlations and the orbital Kondo effect, Sr2(Ba2)V04 under 
pressure - a potential superconductor (477), NiO (32; 478; 479; 481; 480), 
YTiOa (310), TiOCI (482; 483), TlaMnaOy (484), MnO, FeO and CoO (481), 
as well as the Verwey transition in Fe304 (485). 



6.3.1 Ferro-orbital order in LaTiOa- LDA+DMFT calculations started 

with the investigation of Lai-a^Si'x-TiOa by Anisimov et al. (31), who used 
the IPT method as a DMFT solver. Subsequently, Zolfl et al. (279) repeated 
these calculations employing LDA+DMFT(NCA), and Nekrasov et al. (159) 
using LDA+DMFT (QMC). The latter authors also compared the application 
of different DMFT solvers, approximative and numerically exact ones, for a 
realistic material calculation, i.e. doped LaTiOs which is a strongly correlated 
metal close to a Mott-Hubbard metal-insulator transition. The results showed 
that the method for solving the DMFT equation matters, as has already been 
discussed in Section 5.6 with Fig. 19 showing the spectrum of lightly doped 
LaTiOs. Very recently, Craco et al. (312) performed new LDA+DMFT(IPT) 
calculations for LaTiOa, including symmetry breaking. Pavarini et al. (310) 
and Craco et al. (312; 363) reported a ferro-orbital order. This is of partic- 
ular interest since it rules out the orbital-liquid picture of Khaliullin et al. 
(364; 365). The LDA+DMFT results hence provided for an important piece 
of information concerning the controversial debate on the physics of LaTiOs. 

Recently, also SrTiOs/LaTiOs heterostructures have been studied by 
LDA+DMFT and PES (366). Such transition metal oxide heterostructures 
promise to be a vivid area of research in the future, and LDA+DMFT can 
provide for the necessary theoretical support. 



6.3.2 Mott-Hubbard transition in V2O3. A particularly important system 
is V2O3 which undergoes the famous Mott-Hubbard metal-insulator transition 
(367; 368), sec the phase diagram Fig. 26. Held et al. (353; 369) investigated 
paramagnetic V2O3 and Cr-doped V2O3 by LDA+DMFT(QMC), describing 
a Mott transition at a reasonable strength of the Coulomb interaction. The au- 
thors reported reasonable agreement with photocmission spectroscopy (PES) 
experiments by Schramme et al. (370), see Fig. 27, as well as with the ex- 
perimentally expected spin and orbital configuration (371). Later, Mo et al. 
(372; 373) reinvestigated the PES spectrum of V2O3 by new bulk-sensitive 
PES at the Spring-8 synchrotron. Their results represent a big step forward 
on the experimental side since the V2O3 quasiparticle peak could be resolved 
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Figure 26. Left: Phase diagram of V2O3 as a function of Cr and Ti doping and/or pressure P. 
Within the paramagnetic phase a Mott-Hubbard transition between metallic (PM) and insulating 
phase (PI) occurs upon increasing U/W (reproduced from (368)). Right: LDA+DMFT result for 
the effective chemical potential — ReS(a;o) for the a\g and eg orbitals. (reproduced from (369)) 



for the first time. The PES results are also in better qualitative agreement 
with LDA+DMFT(QMC), but still show a broader quasiparticle peak with 
more spectral weight than theoretically expected. 

The LDA+DMFT(QMC) results (353; 369) also include important differ- 
ences to the Mott-Hubbard transition in a one-band Hubbard model. While 
the latter is characterised by the divergence of the effective mass (or vanishing 
quasiparticle weight Z — > 0), the effective mass for the aig band remains fi- 
nite at the transition in V2O3, only that of the eg band diverges. This reflects 
the more complicated nature of the Mott-Hubbard transition in multi-orbital 
systems. As a detailed analysis by Keller et al. (369; 374) revealed, the aig 
orbital becomes insulating since the effective chemical potential ji — ReIl(a;o) 
moves out of the band edges of the non- interacting DOS, see Fig. 26 (right 
panel). Hence close to the Fermi energy, the aig band behaves as if it becomes 
band-insulating at the transition. Shifts of the aig effective chemical potential 
were also reported for BaVSs (375), albeit not leading to insulating aig bands 
in the paramagnetic phase. For BaVSs, these shifts and corresponding changes 
of the orbital occupation explain the the nature of the charge density wave. 

Another important LDA-I-DMFT prediction was that the Mott gap is filled 
with spectral weight upon increasing temperature. This characteristic feature 
has been confirmed recently in PES experiments by Mo et al. (377). It also re- 
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Figure 27. Comparison of the LDA+DMFT(QMC) spectra with photoemission (left; (370; 372)) 
and X-ray absorption (right; (376)). The photoemission experiments are for two different photon 

energies (reproduced from (369)). 
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Figure 28. Same comparison as on the left hand side of Fig. 27, but using the NMTO t2g 
Hamiltonian instead of the LDA DOS and slightly different Coulomb interaction parameters 

(reproduced from (382)). 



fleets in the do eonductivity which was reinvestigated by Limelette et al. (378) 
who found critical exponents which agree with those of the DMFT Mott transi- 
tion scenario (379). Laad et al. (380; 381) also studied the Mott-Hubbard tran- 
sition in V2O3, using LDA-I-DMFT(IPT) and arriving at the conclusion that 
the Mott-Hubbard transition is driven by changes of the trigonal distortion 
under pressure. Using the NMTO-downfolded Hamiltonian instead the simpli- 
fication to the LDA DOS employed in the earlier studies (353; 369; 380; 381), 
Poteryaev et al. (382) recently recalculated the LDA-I-DMFT spectrum and 
reported somewhat better agreement with experiment (377; 383), see Fig. 28. 

In model calculations corrections due to non-local correlations have been 
discussed (384; 138). In particular close to the antiferromagnetic phase transi- 
tion, paramagnon-like excitations result in strong non-local correlations (138). 
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Figure 29. Left: Comparison of single-site and two-site LDA+DMFT spectra for VO2 for the 
(insulating) monoclinic crystal structure, (reproduced from (122)). Right: Comparison of these 
LDA-|-cluster DMFT spectra with photoemission experiments for the rutile (R) and monoclinic 

(Ml) phase (reproduced from (389)). 

Besides, non-local correlations, also the electron-phonon coupling is expected 
to be of relevance close to the Mott-Hubbard transition (385; 386). The strong 
changes of the electronic degrees of freedom couple to the lattice. Hence, 
while the local electronic correlations described by LDA+DMFT are the driv- 
ing force for the transition, we have to expect corrections to the present 
LDA+DMFT results in the immediate vicinity of the Mott-Hubbard tran- 
sition and close to the onset of magnetic order. 



6.3.3 Peierls transition in VO-z- There is a metal-insulator transition in 
another vanadate: VO2 whose high-temperature rutile phase is metallic while 
the low-temperature monoclinic phase is insulating. One can make use of this 
effect in "smart" windows, which become reflective (metallic) if bright sunlight 
heats them up. 

Already within LDA (387), the monoclinic phase is almost gapped due to the 
dimerisation of V atoms. Using single-site LDA+DMFT, VO2 and has been 
studied by Laad et al. (388) and by Liebsch et al. A more realistic scenario 
for the insulating nature of VO2 has been proposed however by Biermann et 
al. (122), using a two site cluster DMFT. These two sites form a spin singlet, 
triggering the insulating behaviour (122). Hence, the gap in VO2 shown in Fig. 
29, is not a Mott-Hubbard gap but a Peierls gap. 
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6.3.4 Orbital selective Mott-Hubbard transition in Co/^—xSr^RuO/^. 
Among the first materials studied by LDA+DMFT was Ca2-a;Sra;Ru04 an 
unconventional superconductor. The LDA+DMFT calculations by Liebsch 
and Lichtenstcin (390) and by Anisimov et al. (391; 392) arc very different 
concerning the size of the Coulomb interaction of the 4d Ru orbitals, i.e. 
U = 0.8 eV, J = 0.2 eV in (390) and U = 1.7 eV, J = 0.7 eV in (392). Both 
groups did not address the superconducting phase, and restricted themselves 
to the normal phase. Actually, the study of superconducting transition metal 
oxides by DMFT is very difficult, since low temperatures are required. Fur- 
thermore, non-s-wave superconductivity is only possible with more formidable 
cluster DMFT calculations, which for the one-band Hubbard model indeed 
show d-wave superconductivity, see Ref. (128) and Ref. (117). These one- 
band calculations have some relevance for the cuprates but will not be dis- 
cussed here as they arc rather model than material calculations. Moreover 
these calculations have already been reviewed in Ref. (115). Coming back 
to Sr2Ru04, Liebsch and Lichtenstcin (390) reproduced the experimental It- 
resolved spectrum very well, including the Hubbard side structure and the 
Fermi surface. Anisimov et al. (391) reported that the Mott-Hubbard transi- 
tion in Ca2-a;Sr3;Ru04 occurs subsequently for different orbitals so that there 
is a region where some orbitals are insulating and others are metallic, coining 
the name orbital-selective Mott transition. This orbital-selective Mott tran- 
sition arises due to different widths of the non-interacting bandwidth of the 
two-different types of t2g bands and has been a subject of intensive model 
studies thereafter (393; 394; 395; 396; 397; 398; 399; 400; 401; 402; 403; 404; 
405; 406; 407; 408; 409; 410), leading to, at first glance, puzzling discrepancies. 
These discrepancies were settled when it was realised that the Mott-Hubbard 
transition is very different for a system with Z2 symmetric Hund's exchange 
(393; 394; 396; 404; 405; 406; 409) and SU(2) symmetric Hund's exchange 
(395; 397; 398; 399; 400; 401; 402; 403; 350; 411; 412). In the former case the 
Kondoesque quasiparticle peak cannot occur in the immediate vicinity of the 
Mott-Hubbard transition since the Z2 symmetric Hund's exchange forms a 
Sz = ±1 spin and a spin-flip from Sz = +1 to Sz = —I is not possible. Con- 
sequently there is a dip in the spectral function for both the narrow and wide 
band, see Fig. 30. Sakai et al. (292) reinvestigated Sr2Ru04 using a SU(2) 
symmetric Hund's coupling. But, for the rather small Coulomb interaction 
of (390) (for which both bands have pronounced quasiparticle peaks) Sakai 
et al. did not yet find pronounced differences between SU(2) and Z2 Hund's 
exchange. 



6.3.5 'Kinks' in SrVOs. Among transition metal oxides, SrVOs is particu- 
larly simple because it has (i) a 3d^ electronic configuration and (ii) a perfectly 
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Figure 30. Orbital selective Mott-Hubbard transition for a two band model system with 
bandwidths 4 and 2 (a.u.) for the two bands respectively. The left side shows the situation for 
SU(2) symmetric Hund's exchange. At U=2.6 (a.u.), the wide band is still metallic (upper left 
panel) whereas the narrow band (lower left panel) is already insulating (reproduced from (403)). 
this is in contrast to the right side, showing the Z2 symmetric situation. Here, the intermediate 
phase has a gap for the narrow band (full line) and a dip in the spectrum at the Fermi energy for 
the wide band (dashed line) (reproduced from (396)). 



cubic perovskite lattice structure, see Fig. 31. Due to (i), the effect of Hund's 
exchange interaction on the ground state properties is less crucial since this 
exchange interaction only takes effect for two or more electrons. The cubic 
symmetry (ii) on the other hand results in three degenerate t2g bands at the 
Fermi energy and allows for the simplification described in Section 4.1.5, i.e. 
using the DOS instead of the full LDA bandstructure. Hence, SrVOa can be 
considered as a transition metal oxide prototype. Despite this simplicity, there 
has been some debate concerning the series Ca^Sri-j^VOs in which Ca doping 
X leaves the Sd^ configuration unchanged but results in an orthorhombic dis- 
tortion. Interest in this 3(i^ series was initiated by Fujimori et al. (413) who 
reported a pronounced lower Hubbard band in PES. While thermodynamic 
properties such as the Sommerfeld coefficient, resistivity and paramagnetic 
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Figure 31. Left: cubic perovskite crystal structure of SrVOs; right: LDA bandstructurc calculated 
using LMTO and NMTO (reproduced from (419; 420)) 



susceptibility were reported to be essentially independent of x (414; 415; 416), 
PES (413) and Bremsstrahlungs isochromat spectra (BIS) (417) suggested 
dramatic differences between CaVOs and SrVOs, leading to the suggestion of 
a Mott-Hubbard transition with increasing Ca doping x (418). 

This puzzling discrepancy was settled through bulk-sensitive (high photon 
energy) PES by Maiti et al. (421) and by Sekiyama et al. (422; 155) which, in 
mutual agreement with LDA+DMFT calculations (423; 155; 310; 419; 420), 
showed similar spectra for CaVOs and SrVOs. Hence the earlier reported 
differences were attributed to the surface-sensitivity of low photon energy 
PES. This is also supported by LDA-I-DMFT calculations of Liebsch et al. 
(424; 425) which show pronounced differences between SrVOs surface and 
bulk spectra. Fig. 31 shows the LDA bandstructurc, and Fig. 32 the calcu- 
lated LDA-I-DMFT(QMC) spectra in comparison with experiment. Note that 
the LDA-I-DMFT spectra are parameter free since the inter-orbital Coulomb 
interaction U = 3.55 eV and Hund'sxchange J = 1.0 eV have been obtained 
through constrained LDA calculations (423; 155). Pavarini et al. (310) also 
systematically studied the similar systems LaTiOs and YTiOs which have 
an increasingly smaller LDA bandwidth so that these materials become indeed 
Mott-Hubbard insulators, see Fig. 33. 

Despite SrVOs being a simple material, the LDA-I-DMFT(QMC) results 
led nonetheless to a surprise: The dispersion of the quasiparticle peak shows 
'kinks' at — 0.25eV and (less pronounced) at 0.25 eV. This was unexpected 
since such 'kinks' are usually associated with phonon modes, particularly for 
high temperature superconductors (426), or the coupling to other bosonic 
degrees of freedom such as spin fluctuations. The LDA-I-DMFT calculations 
showed that such 'kinks' emerge naturally in strongly correlated systems with 
narrow quasiparticle peaks. Let us note that kinks have also been observed 
experimentally for SrVOs by Yoshida et al. in angular-resolved PES (427). 
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Figure 32. Parameter-free LDA-fDMFT(QMC) spectra for SrVOs and CaVOa compared to PES 
and XAS experiments below and above the Fermi energy respectively (reproduced from (419)). 
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Figure 33. LDA-t-DMFT(QMC) spectra (thick Unes) obtained by Ref. (310) for the 3d^ series 
SrVOs, CaVOs, LaTiOs and YTiOs. The thin lines represent the corresponding LDA DOSes 

(reproduced from Ref. (310)). 

Using the DMFT self-consistency equation, Byczuk et al. (428) were able to 
show mathematically how a central peak in the spectrum necessarily results 
in a wiggle in the self energy and, hence, a 'kink' in the dispersion relation. 
Landau's Fermi liquid theory is restricted to the excitations between Fermi 
energy and 'kink'. The rest of the central peak is strongly damped and follows 
another effective mass renormaUsation. That is. 



-E^k = ZpL^k for \E]i\<uj^ 
Ek = Zcp€k±c for|£^k| >^^*, 



(152) 
(153) 



February 2, 2008 



18:28 



Advances in Physics 



ESCwDMFT'ie 



Electronic structure calculations using dynamical mean field theory 



95 




NMTO downfolded V3d{t,^) bands 

Renormalizedby m*/m NMTO bands 

• Quasiparticle bands obtained from LDA+DMFT(QMC) 



G(0,0,0) 



M(lU,7t,0) 



X(7U,0.0) 



G(0,0,()) 



R{7C,7t,7t) 



Figure 34. Comparison of the LDA bandstructurc for SrVOs (full lines) with that from 
LDA+DMFT(QMC) (dots). Basically the LDA+DMFT(QMC) describes renormalised 
quasiparticle bands with weight = 1.9. However, there are deviations in the form of 'kinks' at 
±0.25 eV (the Fermi energy is set to zero; reproduced from Ref. (420)) 

where u;^, is the "kink" energy, the dispersion relation of the interacting 
system, ZpL the Fermi hquid renormahsation factor, Zcp the renormahsation 
factor for the dispersion beyond the Fermi-hquid regime, c is a constant. From 
the bare (LDA) dispersion and one of the renormahsation factors ZpL, Zqp 
and w^K, the other two can be calculated in a simple way (428). This allows 
for example to calculate from the linear specific heat coefficient (~ 1/Zfl) the 
'position of the "kink" and the overall bandwidth in angular resolved PES. 

6.3.6 e'g hole pockets in NaxCO-z- Unconventional superconductivity was 
also discovered in hydrated Naa;C02 (429), resulting in substantial interest 
in this material. This material is difficult to characterise, particularly, the 
0-Co-O angle is affected by the hybridisation and not known. Presently, 
it is also unclear which orbitals are responsible for the superconductivity, 
and even from which orbitals the Fermi surface is composed of. In LDA 
(430; 431; 432; 433; 434; 435) both aig and e'g orbitals have similar centres 
of gravity, leading to the possibility of so-called hole pockets in the e'g bands. 
Also the aig bands seem to play an important role, at least these orbitals 
change most dramatically in LDA when going from normal-conducting mono- 
layer- hydrated cobaltate to superconducting bi-layer-hydrated cobaltate (434) . 
Since these t2g bands are very narrow. Coulomb correlations are expected to 
play a major role. Within LDA-f-U, Zhang et al. (436; 435) reported the dis- 
appearance of the e'g hole pockets. In contrast, more accurate LDA+DMFT 
calculations by Ishida, Perroni and Liebsch (437; 438) predict the contrary, 
i.e. an enhancement of these pockets. These results have been challenged re- 
cently by Marianetti et al. (440), reporting the opposite behavior, i.e., the 
suppression of the e'g hole pockets in LDA-I-DMFT. This would agree with 
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experiments where these pockets seem to be absent. Such conflicting results 
are possible since, in the case of cobaltates, the results are very sensitive to 
small changes of the LDA crystal-field splitting, see Ref. (435) for a discussion. 
Finally, let us also mention two DMFT model studies, concentrating on bad 
metallic behaviour (439) and the Na-induced potential (441), respectively. 



6.3.7 Heavy-Fermion behaviour in LiV204,. Another interesting transi- 
tion metal oxide is LiV204, the first d system signalling heavy Fermion be- 
haviour by, for example, an unusually large specific heat coefficient, one or- 
der of magnitude larger than in other transition metal oxides. Nckrasov et 
al. (442) did realistic LDA+DMFT calculations for LiV204 and investigated 
whether the scenario of Anisimov et al. (443) for the heavy Fermion behaviour 
holds. The basic idea of Ref. (443) is a separation of the partially filled i2g 
electrons into localised ones forming local moments and delocalised ones pro- 
ducing a partially filled metallic band. Then, the hybridisation between those 
two subsets of electrons, as in /-electron materials, can give rise to heavy 
Fermion effects. Nekrasov et al. (442) found a strong competition between 
antiferromagnetic direct (the exchange constant corresponds to an energy ~- 
450 ks K) and ferromagnetic double exchange (« 1090 K). With these esti- 
mates it appears to be reasonable that these two contributions almost cancel 
so that the Kondo exchange (~-630 /cs K) prevails, resulting in heavy Fermion 
Kondo physics. But, since the energy differences are rather small and the 
calculations were performed at relatively high temperatures and with a nu- 
merical error, the results of (442) do not allow for a final conclusion. In their 
LDA+DMFT(IPT) calculations, Laad et al. (444) on the other hand stressed 
the importance of geometrical frustration for the heavy Fermion behaviour 
which forms the basis of another explanations of the unusual behaviour in 
LiV204 (445; 446; 447; 448; 449; 450; 451; 452). 

More recently, Arita et al. (453) managed to investigate the low tempera- 
ture behavior of LiVOa using PQMC as an impurity solver. These calculations 
indeed show a sharp peak above the Fermi energy, in agreement with experi- 
ments (454), see Fig. 35. The physical origin of this sharp peak is however very 
different from the aforementioned scenarios. The sharp peak emerges from the 
physics of the aig orbital, which is a lightly doped Mott-Hubbard insulator, 
i.e., a metal with a large mass renormalisation. 



6.3.8 Colossal magnetoresistance in manganites. Because of the colos- 
sal magnctorcsistance (CMR) (455; 456), manganites such as Lai_3;Caa;Mn03 
have been at the focus of reasearch during the last years, including many 
DMFT studies (74; 58; 110; 111; 457; 458; 75; 459; 460; 112; 461; 462; 463; 464). 



February 2, 2008 18:28 Advances in Physics ESCwDMFT'16 



Electronic structure calculations using dynamical mean field theory 



97 





" 

reference (Au) 


1'* 


" 17K 

37K 

. 50K 





40 20 Ef -20 
Binding Energy (meV) 



40 20 & 

Binding Energy (meV) 



-20 



Figure 35. Left: Spectral function of LiVOs as calculated by LDA+DMFT for the aig and the two 
Cg orbitals at zero temperature (PQMC) and T = = l/40eV (reproduced from (453)). Right: 
Corresponding experimental PES spectrum in the vicinity of the Fermi energy (reproduced from 

(454)). 



In the parent compount (x = 0), Mn has a configuration with cubic symme- 
try so that 3 electrons occupy the lower- lying t2g orbitals, forming a spin 3/2. 
This spin and its coupling to the remaining itinerant electron in the doubly de- 
generate Cg bands constitute the (ferromagnetic) Kondo lattice model, which 
gives rise to ferromagnetism due to the so-called double exchange mechanism: 
A ferromagnetic allignment of the t2g spins is favorable since it maximises 
the kinetic energy of the Cg electrons. This model was solved exactly within 
DMFT by Furukawa (74; 75). Also employing DMFT, Millis and coworkers 
(58; 110; 111) pointed out however that the double exchange of the Kondo 
lattice model is not enough for describing manganites and the CMR in partic- 
ular. The importance of the electron-phonon coupling to Jahn- Teller phonons 
was stressed, and the Kondo lattice model plus Jahn- Teller phonons was stud- 
ied using DMFT (110; 111; 458; 112). This model gives large magnetoresis- 
tances, however only without Ca doping x and also fails to describe other 
experimental aspects. Held and Vollhardt (459) on the other hand stressed 
the importance of electronic correlations induced by the Coulomb interaction 
between the Cg electrons. The arising question, whether the electron-phonon 
coupling or the Coulomb interaction is responsible for the insulating nature of 
the parent compound, was settled recently by LDA-I-DMFT calculations with 
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Figure 36. Optical conductivity cf(u)). Left: theory (reproduced from (464)), right: experiment 
(reproduced from (470)) for the paramagnetic (PM) phase of Lai-^^Sr^MnOa, {n = 1 — x 
electrons/site in the eg orbitals). The dotted Une is without electron-phonon coupUng, showing a 
metaUic Drude peak; the dashed Hne is the optical conductivity for the ferromagnetic phase (FM) 
at a; = 0.2. Inset: The PM resistivity strongly increases with decreasing T, resulting in a 'colossal' 
magnetoresistance when going from the PM to the FM phase. 

a static Jahn- Teller distortion (463), showing that both Jahn- Teller distortion 
and Coulomb interaction are needed: The Jahn- Teller distortion gives rise to 
a crystal field splitting of the two Cg orbitals which is largely enhanced by the 
Coulomb interaction so that a gap emerges. Since the Cg electrons are also 
spin-alligned to the t2g spin, LaMnOa has one spin-polarised Cg electron per 
site and is hence an insulator, albeit not a Mott-Hubbard insulator. With the 
crystal-field splitting being reduced under pressure, Yamasaki et al. (463) were 
also able to explain the pressure-induced insulator-to-metal transition found 
experimentally by Loa at al. (465). Yang and Held (464) subsequently stud- 
ied doped manganites with dynamical Jahn- Teller phonons instead of a static 
distortion. The authors reported the trapping of Cg electrons as Jahn- Teller 
polarons, an effect already described earlier without eg-Cg Coulomb interaction 
(110; 111; 458; 112). However, the electronic correlations strongly enhance the 
tendencies towards polaron formation. The trapping of Cg electrons as polarons 
results in an insulating-like paramagnetic phase with a pseudogap structure in 
the spectrum at the Fermi energy, irrespectively of doping x, in agreement with 
experiment (466; 467; 468; 469). Including the effects of electronic correlations, 
also the optical conductivity (470; 471; 472; 473; 474) could be described, see 
Fig. 36, and a large CMR were reported for the doped system, see the inset 
of Fig. 36. 



6.4 Other Materials 
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Figure 37. Left: Density of states for NiMnSb for spin-up (upper panel) and down (lower pannel). 

Due to electronic correlations, nonquasiparticle states appear within the LSDA gap when using 
DMFT. Right: Spin-polarisation P and magnetization Af as a function of temperature for FeMnSb. 
The nonquasiparticle states lead to a reduction of the spin polarisation (reproduced form (491) and 

(493), respectively). 



6.4.1 Half-metallic ferromagnetism in Heussler alloys. For the tech- 
nical reahzation of spintronic (486; 487) devices, materials with high spin- 
polarisation are necessary. A promising candidate are half-metallic ferromag- 
nets which show metallic behavior for one-spin species and insulating behavior 
for the other one, allowing in principle for a fully spin-polarized current. Such 
materials have been intensively studied by LDA. However, many of the po- 
tential candidates have d electrons so that electronic correlations cannot be 
neglected. An example of such correlation effects are nonquasiparticle states 
in the minority spin-band above the Fermi energy (488; 489; 490). This effect 
was also studied in realistic LDA-I-DMFT calculations for the Heussler alloys 
NiMnSb (491; 492) and FeMnSb (493), see Fig. 37. Due to these nonquasiparti- 
cle states, the polarisation is reduced at finite temperatures if electronic corre- 
lations are strong, see left part of Fig. 37. Since the emerging non- Fermi-liquid 
physics is related to magnons, which can only be described very rudimentary 
by the local DMFT, the inclusion of non-local correlations is necessary for a 
more accurate description. The nonquasiparticle states lead to a reduced po- 
larisationwhich is a severe disadvantage for spintronic sources which optimally 
would be fully polarised. Other candidates for half-metallic ferromagnetism 
which have been studied by LDA-^DMFT are Gai_a,Mna,.As (494; 495), VAs 
(496), CrAs (497) and CrOs (498; 499). 



6.4.2 Superconductivity in A^Cqq. Synthesized fullerenes display a varity 
of unusual properties, including superconductivity up to relatively high tem- 
peratures, e.g. 33K for CsajRbyCgo (500). Such superconducting alkali-doped 
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Figure 38. Tc as a function of electron doping in Ceo fuUerenes for different values of the ratio 
Coulomb interaction over band-width U/W (phonon frequency ojph = 0.24W; electron-phonon 
coupling strength A = 0.6 (reproduced from (114)). 



fullerenes A^^Ceo are molecular solids with a low bandwidth so that both the 
phonon frequencies and the Coulomb interaction are compareable to the band- 
width, and the dimensionless electron-phonon coupling is of order 1. Hence, 
it is important to treat the electron-phonon coupling and the Coulomb inter- 
action on an equal footing which is possible through DMFT. Such realistic 
DMFT calculations have been carried out by Capone et al. (113; 501) and by 
Han et al. (114). These authors found that the superconductivity mediated by 
electron-phonon coupling is surprisingly resistant to the Coulomb interaction. 
Together with the strong electron-phonon coupling this hence explains the high 
superconducting transition temperatures, without the need for a new elec- 
tronic mechanism for superconducitivity. Furthermore, Coulomb interaction 
and Jahn- Teller coupling tend to localise the electrons at the fillings n = 2,4, 
which explains why fillings close to n = 3 are most favorable for supercon- 
ductivity, see Fig. 38. This correlation effect is beyond Eliashberg theory and 
explains the experimental oberserved (502) change of Tc with doping. 



6.4.3 Mott-insulating zeolites. Zeolites loaded with alkali metals such as 
K„Ali2Sii2048 are at first glance an unlikely candidate for strong electronic 
correlations. These materials however form 'superatoms', see Fig. 39, which 
act like hydrogen atoms with well difined s and p orbitals. The superatoms in 
turn crystalise, see central panel of Fig. 39, with a controllable tunneling from 
superatom to superatom. Since this tunneling is rather weak, these systems 
then indeed belong to the class of strongly correlated materials, despite the 
ingredients being K, Al, Si and O atoms. As Arita et al. (503) showed by means 
of LDA and DMFT calculations, K„Ali2Sii2048 is Mott-Hubbard insulating 
even when the nominal doping level (averaged n) is fractional, which explains 
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Figure 39. Left: 'Superatom' (a cage) of the undoped zeoUte (dark blue: Si, Hght blue: Al, dark 
green: oxygen, orange and red: K); middle: crystallisation of the cages; right: DMFT spectrum for 
realistic (unisotropic) hopping parameters showing Mott-insulating behavior for realistic U values. 

(reproduced from (503)) 

the experimental fact that K-doped zeohtes are insulators (504; 505). Zeolites 
and fullerenes show that also materials with s and p valence electrons can be 
strongly correlated because of a reduced hopping. Tehn interesting many-body 
physics emerges. 



7 Summary and outlook 



The physics of materials with strong electronic correlations such as transi- 
tion metal oxides and heavy Fermion systems is characterised by renormalised 
quasiparticles or Mott-insulating behaviour. Conventional electronic structure 
calculations, for example, in the LDA or GW approximation, cannot capture 
this kind of physics and the corresponding energy scales. This became possi- 
ble by incorporating DMFT into realistic material calculations, merging the 
strength of LDA (or GW) as an ab initio approach with that of DMFT to 
deal with local electronic correlations between d or f electrons. LDA+DMFT 
proved to be a breakthrough and is by now a standard approach, albeit stan- 
dardised program packages are still imperative and self-consistent calculations 
where the LDA bands are renormalised due to the DMFT correlations are 
still rare. Also the calculation of physical quantities other than the electronic 
spectrum, the quantity inherently obtained as the Weiss field in DMFT, is 
not yet standard. Since in principle all physical quantities can be obtained, 
there are many opportunities for future work. One should be cautious however 
that non-local correlations, which are neglected in DMFT, are not particularly 
important for the physical quantity at hand. 
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So far, LDA+DMFT has been successfully applied to many materials, with 
the focus naturally on those with strong electronic correlations. Among others, 
the volume, the electronic and the phonon spectrum of Pu, the changes of phys- 
ical properties in Ce through the volume collapse, and the ferromagnetism of Fe 
and Ni, including the famous 6 eV in Ni, have been calculated. Most work was 
devoted to transition metal oxides. Here, the Mott-Hubbard metal-insulator 
transition in V2O3 was investigated and the filling of the Mott-Hubbard gap 
with increasing temperature predicted. The gap in VO2 was identified as a 
spin-Peirls gap, and the origin of the colossal magnetoresistance in manganites 
understood: the paramagnetic phase is insulating-like because of the localisa- 
tion of electrons as lattice polaron, assisted by strong electronic correlations. 
Even the, at first glance, simple system SrVOs was found to show interesting 
physics, i.e. 'kinks'. In contrast to the calculations for transition metals and 
/ electron systems, those for transition metal oxides were restricted to the d 
orbitals crossing the Fermi energy, except for one recent calculation on NiO 
(480). 

A severe drawback of most calculations was that the Coulomb interaction 
has been used as a free parameter or was chosen ad hoc. This is against 
the spirit of ab-initio calculations which, as a matter of principle, should be 
parameter-free. Encouraging in this respect are the calculations for the pro- 
totypical 4/ and Sd materials Ce and SrVOa which show that parameter-free 
LDA+DMFT calculations are possible if the constrained LDA method is em- 
ployed for determining the Coulomb interaction. In the case of GW+DMFT, 
the Coulomb interaction is even inherently obtained from GW. Here, the chal- 
lenge is to do DMFT calculations with a frequency-dependent (screened) in- 
teraction. 

One should always keep in mind that LDA+DMFT is an approximation for 
solving the solid state Hamiltonian. The starting point, LDA, is already an ap- 
proximation even for the extended s and p orbitals. The second approximation 
involved is the selection of the 'interacting' orbitals which by principle is basic 
dependent. Because of this also the determination of the Coulomb interaction 
has a considerable uncertainty, less so for well localized 4/ orbitals but more 
for 2>d and particularly Ad orbitals. Similarly, the underlying spectral density 
functional theory is basis dependent, in very contrast to standard density func- 
tional theory. Thirdly, DMFT is a many-body approximation for dealing with 
strongly correlated electron systems. It neglects non-local correlations. More- 
over, the solution of the DMFT impurity problem involves possibly a fourth ap- 
proximation if for example the non-crossing or iterated perturbation theory is 
employed. While numerically exact quantum Monte Carlo (QMC) simulations 
are also possible, they are restricted to temperatures at or above (roughly) 
room temperature since they become too expensive in terms of CPU time at 
lower temperatures. In this respect density matrix renormalization group cal- 
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culations, continuous-time and projective quantum Monte Carlo simulations 
are promising alternatives for zero (or low) temperature. 

A major challenge for the future will be to go beyond DMFT, i.e. taking 
into account non-local correlations on top of the local DMFT correlations in 
electronic structure calculations. In this respect, cluster extensions of DMFT 
have already been applied intensively to the two-dimensional Hubbard model 
with a focus on unconventional superconductivity. The numerical effort re- 
stricts these cluster extensions to short range correlations within a (relatively 
small) cluster. Realistic multi-band calculations were hitherto even restricted 
to only two sites. This puts particular importance onto the correlation within 
this pair of sites and is hence naturally only appropriate for particular sys- 
tems such as VO2 where such pair correlations dominate. Many important 
phenomena in correlated electron systems on the other hand stem from long 
range correlations, e.g. magnons, the interplay between antiferromagnetic fluc- 
tuations and unconventional superconductivity, quantum critical points, and 
critical behavior in general. How these phenomena arise is quite well under- 
stood diagrammatically for weak coupling, which as a matter of course is not 
the proper starting point for strongly correlated electrons. But if we replace 
the bare Coulomb interaction of the weak coupling perturbation theory by the 
local (fully irreducible) vertex, we generate a set of diagrams which includes 
both the local DMFT diagrams and the dominant diagrams for long-range cor- 
relations. This is the very idea of the dynamical vertex approximation which 
has been introduced most recently along with alternative approaches. An- 
other long term goal for electronic structure calculations with DMFT is the 
calculation of the the ionic positions hitherto taken from experiment, the in- 
clusion of lattice dynamics, and eventually molecular dynamics on the basis 
of LDA+DMFT. 

Continuing the fruitful cooperation between bandstructure and many-body 

physicists, wc can optimistically face these challenges in the future. Electronic 
structure calculations employing DMFT or its diagrammatic and cluster ex- 
tensions will be a prospering area of research. In particular, such approaches 
will further improve our understanding of systems with strongly correlated 
electrons. The advances in electronic structure calculations through DMFT 
put our ability to predict physical quantities of such strongly correlated ma- 
terials onto a similar level as conventional electronic structure calculations for 
weakly correlated materials — at last. 
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